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We review some recent (mostly ours) results on the Anderson localization of light and electron 
waves in complex disordered systems, including: (i) left-handed metamaterials, (ii) magneto-active 
optical structures, (iii) graphene superlattices, and (iv) nonlinear dielectric media. First, we demon- 
strate that left-handed metamaterials can significantly suppress localization of light and lead to an 
anomalously enhanced transmission. This suppression is essential at the long- wave length limit in 
the case of normal incidence, at specific angles of oblique incidence (Brewster anomaly), and in the 
vicinity of the zero-£ or zero-// frequencies for dispersive metamaterials. Remarkably, in disordered 
samples comprised of alternating normal and left-handed metamaterials, the reciprocal Lyapunov 
exponent and reciprocal transmittance increment can differ from each other. Second, we study 
magneto- active multilayered structures, which exhibit nonreciprocal localization of light depending 
on the direction of propagation and on the polarization. At resonant frequencies or realizations, 
such nonreciprocity results in effectively unidirectional transport of light. Third, we discuss the 
analogy between the wave propagation through multilayered samples with metamaterials and the 
charge transport in graphene, which enables a simple physical explanation of unusual conductive 
properties of disordered graphene superlatices. We predict disorder-induced resonances of the trans- 
mission coefficient at oblique incidence of the Dirac quasiparticles. Finally, we demonstrate that 
an interplay of nonlinearity and disorder in dielectric media can lead to bistability of individual 
localized states excited inside the medium at resonant frequencies. This results in nonreciprocity of 
the wave transmission and unidirectional transport of light. 
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I. INTRODUCTION 

Anderson localization is one of the most fundamental 
phenomena in the physics of disordered systems. Be- 
ing predicted in the seminal papei^ for spin excitations, 
and then extended to electrons and other one-particle 
excitations in solids^ ^ and classical waveJ^K^, it became 
a paradigm of the modern physic^. The study of this 
phenomenon remains a hot topic throughout its more 
than 50-years history. It is constantly stimulated by new 
experimental results, including the most recent observa- 
tions in microwaveJ^HI^, opticJ^^H^l, and Bose-Einstein 
condensates^^. 

Being a universal wave phenomenon, Anderson local- 
ization has natural implications in novel exotic wave 
systems, such as photonic crystals, meta- and magne- 
tooptical materials, graphene superlattices. Indeed, left- 
handed metamaterials, nonlinear and magnetooptical 
materials, and graphene L3M^ are involved in design and 
engineering of various multilayered structures operating 
in a broad spectral range, from optical to microwave fre- 
quencies. Random wave scattering and localization natu- 
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rally appear in such systems, either due to technological 
imperfections or owing to the intentially designed ran- 
dom lattices. Importantly, exotic properties of the con- 
stituent materials essentially require consideration of the 
interplay of the Anderson localizati on w ith various ad- 
ditional effects: absorptio n and ga iii^^^^^^, polarization 
and spin^i^, nonlinearity^^^^^^^^H^, and magnetoopti- 
cal phenomengP^M^. In this review, we describe novel 
remarkable features of Anderson localization of waves 
in multilayered structures composed of non-conventional 
materials with unique intrinsic properties. 

We start our review with Sec. [TTl which introduces the 
basic concepts and general formalism describing the wave 
propagation, scattering, and localization in in random- 
layered media. Anderson localization originates from 
the interference of multiply scattered waves, manifest- 
ing itself most profoundly in one-dimensional (ID) sys- 
tems where all states become localized^^ Due to 
one- dimens ional geometry, such systems are well ana- 
Iy2g(j2j39i40|^ including the mathematical level of rigorous- 
ness of the results^^^^. We describe the exact transfer- 
matrix approach to the wave propagation and scattering 
in layered media. The main spatial scale of localization, 
i.e., localization length, can be defined in two ways: (i) via 
the Lyapunov exponent of the random system and (ii) via 
the decrement of the wave transmission dependent on the 
system. In usual Anderson-localization problems, these 
two localization lengths coincide with each other. 

In Section IIIII we consider transmission and localiza- 
tion properties of the multilayered H-stacks comprised of 
normal materials with right-handed i?— layers and mixed 
M-stacks, including also left-handed L— layers with neg- 
ative refractive inde:jP^. The opposite signs of the phase 
and group velocities in metamaterials lead to partial or 
complete cancellation of the phase accumulation in mul- 
tilayered M-stacks. We show that this cancellation sup- 
presses the interference of multiple scattering waves and 
the localization itsel P^MMI. Using the weak scattering 
approximation (WSA^P^EH^ give detailed analytical 
and numerical description of transmission and localiza- 
tion properties of both M- and H-stacks and reveal a 
number of intriguing results. Namely: (i) in the long 
wave limit localization lengths defined via the Lyapunov 
exponent and transmission decrement differ from each 
other in M-stacks, (ii) there exist two ballistic regimes 
in the H-stacks, (iii) essential suppression of localization 
at special angles in the case of oblique incidence (Brew- 
ster anomaly) and in the vicinity of special frequencies 
(zero-e or ze ro-/i frequencies) is observed. Finally, in 
Section [mG] we discuss an anomalous enhancement^^ of 
wave transmission in minimally disordered alternated M- 
stacks of metamaterials, where the layer thicknesses are 
equal and only dielectric permittivities (or only magnetic 
permeabilities) vary. 

Section [TV| is devoted to the study of novel localization 
features in novel materials. We start with discussion of 
localization of light propagating through magneto-active 
multilayered structures, with either Fara day or Cotton- 
Muton (Voigt) geometries (Section IV A). We show that 
magnetooptical effects can significantly affect the phase 



relations, resulting in nonreciprocal localization depend- 
ing on direction of the wave propagation and polariza- 
tion of light. At resonant frequencies corresponding to 
the excitation of localized states inside the sample, a 
nonreciprocal shift of the the resonance results in effec- 
tively unidirectional transmission of light^^. In Section 
|IVB[ conducting properties of a graphene layer subject to 
stratified electric field are considered. The close analogy 
between charge transport in such system and wave trans- 
mission through multilayered stack^^ underpins remark- 
able conductive properties of disordered graphend^^. We 
predict disorder-induced resonances of the transmission 
coefficient at oblique incidence of electron waves. Finally, 



in Section IV C we examine the interplay between non- 
linearity and disorder in resonant transmission through a 
random- leyered dielectric medium^^. Owing to effective 
energy localization and pumping, even weak Kerr non- 
linearity can play a crucial role leading to bistability of 
Anderson localized states inside the medium. Akin to 
the magneto-optical structures, this brings about unidi- 
rectional transmission of light. 



II. RANDOM MULTILAYERED STRUCTURES 
A. Transmission Length and Lyapunov Exponent 

As it was mentioned above, ID Anderson localization 
results is exponential decay of the transmission coefficient 
with the length L of the sample. For multilayered sys- 
tems, it worth to use the total number of layers N and 
mean layer thickness L/N. In what follows we use di- 
mensionless variables measuring all lengths in mean layer 
thickness units while the time dependence is chosen in the 
form e"*'^^. For simplicity throughout all this review we 
mainly consider the lossless stacks. The detailed results 
concerning to the case of stacks with losses can be found 
in original works. 

Introduce the dimensionless transmission length 1^ on 
a realization 
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and " averaged" TV-dependent dimensionless transmission 
length It = It{N) of a multilayered A/"— layered stack 
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Here T/v is the stack amplitude transmission coefficient 
related to its transmittivity T/v by equality T/v = ITat^. 
Due to self- averaging of In |T/v|/A^, both these lengths It 
and In tend to the same limit 



lim Iq 

N^oo 



lim In 

N^oo 



I, 



(II.2) 



as the number N of layers tends to infinity. FollowinJ^ 
we recall / as localization length. This localization length 
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is related directly to the transmission properties. Its 
reciprocal value is nothing but decrement of the stack 
transmission coefficient. 

Transmission coefficient entering these equations is 
naturally expressed in terms of the total T-matrix of the 
stack written in the running wave basis. Consider trans- 
mission of the plane wave incident normally from the left 
to the stack comprised of even number N of layers and 
embedded into free space. In the simplest case, the wave 
is described in terms of two component vector of, say, an 
electric field e. Within a uniform medium with dielectric 
permittivity s and magnetic permeability /i, the field e 
has the form 



e(z) 



k = -y/siJ, 
c 



(11.3) 



with 2;-axis directed to the right (here and below all 
lengths of the problem are dimensionless and measured 
in the mean layer thickness). 

If the components of vector e are norm alized by such a 



way that the energy flux of the wave (II. 3 ) is |e+ p - 
then the amplitudes 



eL,R 



(11.4) 



of the field from both sides out of the A''— layer stack are 
related by its transfer matrix T{N) 



h 



{\Rn?) 
2N '■ 



(ILi 



involving t he ave rage reflectance^^. This follows directly 
from Eq. (II.l) by virtue of the current conservation 
relationship, \Rn\^ \Tn\^ = 1. The length h in this 
equation is termed the ballistic length. 

Accordingly, in studies of the transport of the classical 
waves in one-dimensional random systems, the following 
spatial scales arise in a natural way: 

• It — the transmission length of a finite sample 
(ILll), 



• / — the localization length (II. 2) related to trans- 
mission properties, and 



• h — ballistic length (II.8). 



The exponential decrease of transmission coefficient 
with the stack size is only manifestation of Anderson lo- 
calization. The phenomenon of localization itself is the 
localized character of eigenstates in infinite disordered 
system with sufficiently fast decaying correlations. The 
quantitative characteristic of such a localization is the 
Lyapunov exponent which is increment of the exponen- 
tial growth of the currentless state with a given val ue at 
certain point far from this point. The amplitude (II. 4) 



of the currentless state in the inhomogeneous medium in 
the basis of running waves can be parameterized as 



e\L = f{N)(^R, 



(11.5) 



which is expressed via transmission and reflection coeffi- 
cients of the stack as 



f (TV) 
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where asterisk stands for the complex conjugation. 
The methods of calculation of transmission coefficient 



T{N) = (fn)' 



(11.7) 



are discussed in the next Subsection. 

In what follows, we consider stacks composed of weak 
scattering layers with reflection coefficients of each layer 
much smaller than 1. In spite of this, for a sufficiently 
long stack the transmission coefficient is exponentially 
small |T/v| ~ exp {—f<iN) with decrement coinciding with 
reciprocal localization length n = l^^ (localized regime). 
However a short stack comprising a comparatively small 
number of layers is almost transparent \Rn\^ ^ 1 (bal- 
listic regime). Here the transmission length takes the 
form 



R 



-iO 



(II.9) 



where R{z) and 0{z) are the modulus and the phase of 
the considere d cu rrentless solution correspondingly. 

It is knowiP^ that at given initial values ^(0) (i?(0)), 
and 6>(0), the function ^(z) at a sufficiently far point is 
approximately proportional to its distance from the ini- 
tial point. In discrete terms, with the probability 1 the 
positive limit exists 
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^= ^i-S' (n.io) 



lim 



which is called Lyapunov exponent. Its reciprocal value 
we also call localization length 
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(11.11) 



however index ^ reminds that this localization length is 
defined through Lyapunov exponent. 

To compare the two localization lengths / and we 
consider first the continuous case were corresponding dy- 
namical variable ^{z) depends on continuous coordinate 
z. In this case, transmittance of the system with length 
L is exactly expressed a4^^ 
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where ^c(^) and ^s(^) are two independent solutions satis- 
fying so called cosine and sine initial conditions Oc{0) = 
and 6^(0) = 7r/2 and having the same limiting behavior 



from right to left, so that all odd layers j = 2n — 1, are 
of type "a" and all even layers j = 2n are of type , 
n = 1,2, ...,A^/2 (see Fig. [T]). In general case the j— th 
layer is characterized by its dimensionless thickness dj, 
dielectric permittivity Sj and magnetic permeability /ij. 



,= l=limM^=lim m 



(11.13) 



Equations (11.12) and ( |II.13| ) evidently show that in con- 
tinuous case / and exactly coincide. 

In the discrete case (multilayered stack), corresponding 
expression for transmittance reads 

2e^^W+^«W sin {Oc{N) - 6>,(7V)))"^ . (11.14) 

Here the last term in denominator differs from that in 
Eq. (11.12). Moreover, it can change its sign and gen- 



erally speaking can essentially reduce the denominator 
itself thus enlarging transmittance and as a result en- 
larging localization length in compare to /. Thus, Eqs. 
(11.14) and ( |II.1Q ) enable us to state only that I > in 
contrast to the continuous case where these two localiza- 
tion lengths always coincide. In spite of that, studying 
of localization in normal disordered multilayered stacks 
did not show any difference in the two lengths. We will 
see below that such a difference really manifests itself in 
the alternated metamaterial stacks. 

In this review we are mainly interested in the transmis- 
sion length It- This quantity can be found directly by 
standard transmission experiments. At the same time, it 
is sensitive to the size of the system and therefore is best 
suited to the description of the transmission properties in 
both the localized and ballistic regimes. More precisely, 
the transmission length coincides either with the localiza- 
tion length / or with the ballistic length 6, respectively in 
the cases of comparatively long stacks (localized regime) 
or comparatively short stacks (ballistic regime). That is, 

/ N'>1 
b N <^b. 



B. Transfer Matrices and Weak Scattering 
Approximation 

In this Subsection we describe some methods used for 
calculation of transmission length and other transmission 
or/and localization characteristics in various regimes. All 
of them are based on various versions of transfer matrix 
approach. 

Consider the M-stack alternatively comprised of even 
number N of uniform layers labeled by index j = 1, A/" 



Papa 



N N-1 



2n 2n-l 



4 3 2 1 



FIG. 1: (Ref. |43]) Two component multilayered alternative 
stack. 



The total transfer matrix (II.6) is factorized to the 
product 



f(Ar) = £^f^_i...£2fi 



(11.15) 



of the layer transfer matrices tj. 

Note that for considered alternated stack, it is natural 
to join each pair of subsequent layers with numbers j = 
2n — l and j = 2n into one effective cell number n. Then 
the total transfer matrix factorize s to t he product of N/2 
transfer matrices of separate cellJ^^^^^^'^. 

Parameterizing the transfer matrix of the j— th layer 
by its transmission tj and reflection rj coefficients of a 
corresponding layer we obtain the recurrence relations 



1 - Rj-iVj 

' ' i?o = 0, 



1 - R-i-iri 



(11.16) 
(11.17) 



where Tj and Rj are transmission and reflection coeffi- 
cients of the reduced stack comprised of only j first lay- 
ers. These relations provide an exact description of the 
system and will be used later for direct numerical simula- 
tions of its transmission properties. Another possible but 
less effective way is related to direct numeric calculation 
of the total transfer matrix (11.15). 



Relations (11.16) and (11.17) serve as a starting point 



for the weak scattering approximation (WSA) elaborated 
iiP^ and based on assumption that the reflection from 
a single layer is small i.e., \rj\ <C 1. This demand is 
definitely satisfied in the case of weak disorder. Within 
WSA, instead of exact relations (11.16), (11.17) we use for 



the transmission length the following first order approx- 
imations 



Rj 



InTi +lnt^- 



Rj—it^ 



J 



2,3,...,iV. 



(11.18) 
(11.19) 



Note that in deriving Eq. (11.19), we omit the first-order 
term R'j_-^t'jrj. This is uncontrolled action. The omitted 
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term contributes only to the second order of In Tj already 
after the first iteration for not very large number of layers 
j. For sufficiently large j, it should be taken into account. 
Nevertheless as we will see below, this approximation is 
excellent in all wavelength region. 

Neglecting the last term in the right hand side of Eq. 
(11.18) we come to the so-called single-scattering approx- 



imation (SSA), which implies that multi-pass reflections 
are neglected so that the total transmission coefficient is 
approximated by the product of the single layer transmis- 
sion coefficients as well as total transmittance is approx- 
imated by the product of the single layer transmittances 
that results in 



N 



In IT. 



N 



^ln|t, 



In the case of very long stacks (i.e., as the length 
N oo), we can replace the arithmetic mean, 

N~'^ XljLi 1^ \ by its ensemble average (In On the 
other hand, in this limit the reciprocal of the transmission 
length coincides with the localization length. Using the 



energy conservation law. 



1, which applies 



in the absence of absorption, the reciprocal localization 
length in single-scattering approximation may be written 
as 



Using the group property of the interface transfer ma- 
trices: F^" = F^Of^", and F^^ = F^^FO/^ ^^^^i 



transfer matrix is factorized to the product (11.15) where 
the layer transfer matrices are 



Such a representation is especially efficient in the short- 
wave limit where the total transmission coefficient re- 
duces to the product of the transmiss ion co efficients of 
only interfaces (see Ref.^^ and Section IV A). 



Come now to application of the transfer matrix ap- 
proach to calculation of the Lyapunov exponent 7- De- 
fine for each layer the curentless vector ej by Eq. (II. 9) 
with the corresponding values 
Lyapunov exponent is written as 



and Oj. 



In this terms 



7 = lim ^ = lim {^j - 



(11.22) 



(we used Shtolz theorem). The vectors Cj and ej_i satisfy 
the equation 



(11.23) 



and is proportional to the mean reflectance of a single 
random layei*^^. 

The version of transfer matrix approach described 
above is based on consideration of a single layer embed- 
ded into vacuum. This version and related WSA were 
used iipUEH for analytical and numerical study of meta- 
material M-stacks (see Section [lll| ) . 

Another version used iiP^ (Section IV A ) is based on a 



separation of wave propagation inside a layer and through 
the interface between layers (see e.g. Ref.^^). Here wave 
propagation inside the j-th layer is described by diagonal 
transfer matrix 



diag(e 



(11.20) 



where ipj = kjdj is the phase accumulated upon the wave 
propagating from left to right through the j-th layer, and 
^3 ~ ^V^JH- "^^^ interfaces are described by unimodu- 
lar transfer matrices ^ ^"/^ ^ ^ ^/30 corresponding, 
respectively, to transitions (all from left to right) from 
vacuum to the medium 'a', from the medium 'a' to the 
medium from the medium to the medium 'a', and 
from th e med ium 'a' to vacuum. Thus, the total transfer 
matrix (II.6) of the structure in Fig. [l] is 



T{N) = F ^FnSnFn-iSn-1 ... F2S2F1S1F' 



aO 



F2n- 



a/3 



n ■■ 



1,2,..., N/2. (11.21) 



Therefore the difference in the r.h.s. of Eq. ( II.22[ ) is 
some function of 6j-i 



(11.24) 



which explicit form is determined by Eq. 11.23 Using the 
self averaging of the ratio ^ /j and the fact that the phase 
Oj stabilizes^, we finally obtain for Lyapunov exponent 



(11.25) 



where average in the r.h.s. is taken over stationary dis- 
tribution of the phase 0. 

Continuous version of this result was obtained in^ (see 
Eq. (10.2)). Its dis crete version in slightly different terms 
(see Section IIIG) was obtained in^ Note that due to 
existence of the closed formula (11.25) for Lyapunov ex- 
ponent, the task of analytical calculation of the localiza- 
tion length = 7~^ is a simpler problem than that of 
transmission length It- 

The next steps are standard (see e.g. Refs. |2|58) ): 



using (11.23) to get the dynamic equation for the phase 



^, write down corresponding Fokker-Planck equation for 
its distribution, solve it and calculate the average (11.25). 



Moreover, in weakly disordered systems, only the first 
and the second order terms should be accounted for in 
the dynamic equations^ For minimally disordered M- 
stacks defined in Section |Tj this pr ogram was successfully 
realized in^^ (see Section [Lll G below). 
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III. SUPPRESSION OF LOCALIZATION IN 
METAMATERIALS 

Over the past decade, the physical properties of meta- 
materials and their possible applications in modern op- 
tics and microelectron ics, have received considerable at- 
tention (see e.g. Ref d^ ' ^^ ' ^"'^^ '). The reasons for such 
an interest are unique physical properties of metama- 
terials including their ability to overcome the diffrac- 
tion limit"^^^, potential role in cloaking^^, suppression of 
spontaneous emission rate^^, the enhancement of quan- 
tum interference^, etc. One of the first study of the 
effect of randomness^ revealed that weak microscopic 
disorder may lead to a substantial suppression of the 
wave propagation through magnetic metamaterials over 
a wide frequency range. Therefore the next problem was 
to study localization properties of disordered metamate- 
rial systems. 

It was known that, in normal multilayered systems 
comprising right-handed media, the localization length 
is proportional to the square of wavelength A in the 
long- wavelength limit, tends to a constant in a short- 
wavelength regim e, and o scillates irregularly in the inter- 
mediate region^ ^ Natural question arises: how 
inclusion of metamaterial layers influences the localiza- 
tion and transmission effects. 

The study of localization in metamaterials was started 
in Ref. [68 where wave transmission through an alter- 
nating sequence of air layers and metamaterial layers of 
random thicknesses was studied. Localized modes within 
the gap were observed and delocalized modes were re- 
vealed despite the one-dimensional nature of the model. 
Then comprehensive study of transmission properties of 
M-stacks was done iipMH^ Here anomalous enhancement 
of the transmission through minimally disordered (see 
Section |l| M-stacks was revealed^^, non-coincidence of 
the two localization lengths / and was established^, 
polarization^^ and dispersion^^ effects in transmission 
were studied. 

Scaling laws of the transmission through a similar 
mixed multilayered structure were investigated in Ref. 
[69 . It was shown that the spectrally averaged trans- 
mission in a frequency range around the fully transpar- 
ent resonant mode decayed with the number of layers 
much more rapidly than in a homogeneous random slab. 
Localization in a disordered multilayered structure com- 
prising alternating random layers of two different left- 
handed materials was considered in Ref. [70]. Within 
the propagation gap, the localization length was shorter 
than the decay length in the underlying periodic struc- 
ture (opposite of that observed in the random structure 
of right-handed layers). 

Detailed investigation of Lyapunov exponent (and 
therefore localization length l^) in various multilayered 
metamaterials was presented in^^. In the weak dis- 
order limit, explicit expressions for Lyapunov exponent 
valid in all region of wavelengths for various kinds of 
correlated disorder were obtained^^ and analytical ex- 
planation of anomalous suppression of localization was 
dond^ 



Dispersion effects in M-stacks comprised by metama- 
terial layers separated by air layers with only positional 
disorder were considered iiP^^^. Here essential suppres- 
sion of localization in the vicinity of the Brewster angle 
and at the very edge of the band gap was revealed^^, in- 
fluence of both quasi-periodicity and structural disorder 
was studied^^ and effects of some types of disorder cor- 
relation on light propagation and Anderson localization 
were investigated''^. 

In this Section we consider suppression of localization 
in sufficiently disordered M-stacks. In the first four Sub- 
sections we consider the model with non-correlated fluc- 
tuating thicknesses and dielectric permittivities. This 
model possesses the main features cause by the presence 
of metamaterials and at the same time remains compar- 
atively simple. The results concerning disorder correla- 
tions can be found in papers mention in the previous 
paragraph and detailed recent survey The presenta- 
tion is mostly based on works^^ 



A. Model 

We start with the model described at the beginning 
of Suection II B and displayed in Figjl] Electromagnetic 
properties of the j-th layer with given dielectric permit- 
tivity Ej and magnetic permeability /j^j , are characterized 
by its impedance Zj and refractive index Uj 



(IH.l) 



Being embedded into vacuum, each layer can be de- 
scribed by its reflection and transmission coefficients with 
respect to wave with dimensionless length A incident from 
the left 



1 - p2e2ift ■ 



(III.2) 



Here pj = {Zj — l)/{Zj + 1) is Fresnel coefficient, 
= kdjVj^ and /c = 27r/A is dimensionless wavenumber. 

Within our model, dielectric permittivity, magnetic 
permeability and thickness of the j-th layer have the 
forms 



^i = {-iy, d, = i + 5,^, (III.3) 

so that corresponding impedance and refractive index are 



V, = (-1)^(1 + <5(^)). 



(111.4) 
(III.5) 



The thickness fluctuations 5^^^ are independent identi- 
cally distributed zero- mean random variables, as well as 
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all refractive index fluctuations d^^\ To justify the weak 
scattering approximation, we assume that all these quan- 
tities 5^^'^^ are small. 

The considered model possesses some symmetry: sta- 
tistical properties of the fluctuations and absorption coef- 
ficient are the same for L and R layers. As a consequence 
of this symmetry, the scattering coefficients of R and L 
layers are complex conjugate tr = and = r^*, that 
results in the relations 



{giU)) = {g{ti)r, (girr)) = {g{n))*. (m.e) 

valid for any real- valued function g in either the lossless 
or absorbing cases. In more general models this symme- 
try can be broken. 

The model with two parameters (here - thickness and 
refractive index) is in a sense the simplest sufficiently 
disordered model. Further simpliflcation where only one 
of these quantities is random qualitatively changes the 
picture. Indeed, the case of M-stack with only thick- 
ness disorder in the absence of absorption is rather triv- 
ial: such stack is completely transparent (a consequence 
of Zj = 1). On the other hand, M-stack with only 
refractive-index disorder as it was revealed in^^, mani- 
fests a dramatic suppression of Anderson localization - 
essential enlightenment in the long wave region. This in- 



triguing case is considered below in Section III G So here 
we focus on the case where both two types of disorder are 
simultaneously present. 

Speciflc features of transmission and localization in the 
M-stacks look more pronounced in comparison with those 
of homogeneous stack (H-stack) comprised of solely ei- 
ther right-handed or left-handed layers. Therefore albeit 
localization in disordered H-stacks wi th right-han ded lay- 
ers has been studied by many author j^l ^^ * * ^^ * ^ ^, we also 
consider this problem here in its most general formula- 
tion. This consideration enables us to compare localiza- 
tion properties of M- and H-stacks. To describe a H-stack 
composed of onl y R (L ) layers, all multipliers ( — 1)-^ in 



Eqs. (IIL3) and (IIL5) should be replaced by 1 (-1). 



B. Mixed Stack 



Within the version (11.18), ( II.19[ ) of weak scattering 
approximation, contributions from the even and odd lay- 
ers are separated. As a result the transmission length of 
a flnite length M-stack may be cast in the form^ 



< In |t| >, < r >, and < > composed of transmission t 
and reflection r coeflicients of a single right-handed layer: 



-{ln\t\) 



|(r)p+Re((r)^(t^)* 



(IH.9) 



2/1 



b I 1 - exp(-2//) 
'|(r)p+Re((r)2(t2)*) \{r)f 



I 



\n\{t^)\' 



(IH.IO) 



These Eqs. (HI. 7) - ( |III.1Q| ) are valid in the presence of 
absorption. However below to make our treatment more 
transparent, we consider the lossless case. 

The characteristic lengths /, 6, and I are functions 
of wavelength A. The flrst two always satisfy the in- 
equality 1{X) > 6(A), while in the long wavelength re- 
gion the crossover length is the shortest of the three, 
b{X) > 1{X). In the case of a flxed wavelength A, for 
comparatively short stacks with N <C 1{X) the function 
f{N^l) ^ 1, while for sufficiently long stacks N ^ 
it tends to zero f{N, I) ^ 0. Correspondingly, transmis- 
sion length coincides with ballistic length It{X) ~ b{X) 
for short stacks N <C 1{X) and with localization length 
It{X) ~ /(A) for long stacks N 1{X) with the tran- 
sition between the two ranges of N being determined 
by the crossover length [(A). Thus ballistic regime oc- 
curs when the stack is much shorter than the crossover 
length N <C /(A). The localized regime is realized for the 
stacks longer than localization length N ^ ^A). For the 
stacks with intermediate sizes 1{X) ^ N < 1{X)^ transmis- 
sion length coincides with localization length, however 
they correspond to the transition region between ballis- 
tic regime and localized one. 

Alternatively we can consider the stack with a given 
size N and use the wavelength as the parameter gov- 
erning the localized and ballistic regimes. To do this, 
we introduce two characteristic wavelengths, Xi{N) and 
A2(A/'), deflned by the relations 



N = liXiiN)), N = l{X2iN)). 



(III.ll) 



where 

fix) = i^^. (III.8) 

X 

Localization length ballistic length 6, and crossover 
length I are completely described by the three averages 



In these terms, the localized regime occurs if A <C Ai(A^), 
while in the long wavelength region, A X2{N), the 
propagation is ballistic. Intermediate range of wave- 
lengths, Xi{N) < X < A2(A/'), corresponds to transition 
region between the two regimes. 

Consider now example of rectangular distribution, 
where the fluctuations Sj^^ and Sj^^ are uniformly dis- 
tributed over the intervals [—Qu^Qu] and [—Qd^Qd] re- 
spectively and have the same order of magnitude Q^, ~ 
Qd so that the dimensionless parameter 



8 



(11.11). Within the considered model (III.2) - (III.5), its 



,01 
'Ql 



is of order of unity. 

At the next step, we calculate the averages < In \t\ >, 



< |r| >, and < |r | > with the help of Eqs. (IIL2) 



(IIL5), substitute them into Eqs., (III. 9) and (Ill.lO)and 



neglect the contribution of terms of order higher than 
Q^. The resulting general expressions for localization, 
ballistic and crossover lengths are rather cumbersome so 
we present here only their asymptotical forms. 

In the short wavelength region, the main contribution 
to localization length is related to the first term in the 
r.h.s. of Eq. (III. 9) corresponding to the single scattering 



approximation and the localization length is 



/(A) 



12 



A<C1. 



(III.12) 



This means that the size N of the short stack NQ^ <C 1 
is always smaller than localization length and the short 
wave transmission through short stack is always ballistic. 

Opposite limiting case NQl ^ 1 corresponds to the 
long stacks. Here both two regimes are realized and 
transition from localized propagation to the ballistic one 
occurs at the long wavelength A ^ QvVN ^ 1. In- 
deed, asymptotical expressions for all three characteristic 
lengths read 



1{X) 



3A2 3 + C 



3A^ 1 

27r2Q2 4(3 + ^)' 



(III. 13) 



(III.14) 



and 



6(A) 



3A2 



(III.15) 



Note that the single scattering approximation for lo- 
calization length fails in the long wave limit because both 
two terms i n the r .h.s. of Eq. (III. 9) contribute to the 
asymptotic ( |III.13 ). 

Thus in the symmetric weak scattering case, ballis- 
tic, localization, and crossover length in the long wave 
region differ only by numerical multipliers, satisfy the 
inequality < b{X) < 1{X) mentioned above, and 

are pro portional to A^. Two characteristic wavelengths 
(III.ll) correspo nding t o localization length (III. 13) and 
crossover length (III. 14), are proportional to Q^7\Z/V, dif- 



fer only by a numerical multiplier and satisfy the in- 
equality Ai(A^) < A2(A^). For the sufficiently long stacks 
NQl ^ 1, they are lying in long wave region A 1^2 ^ 1- 

Localization properties of infinite stack are described 
by Lyapunov exponent (11.10) or by localization length 



long wave asymptotic calculated with the help of well 
known transfer matrix approach reads 



7 = 



2A2 



(1 + J(^))^ (III.16) 



In the case of rectangular distributions of the fluctua- 
tions of the dielectric constants and thicknesses described 
above, reciprocal Lyapunov exponent reduces to 



/^(A)=7-^(A) 



3A2 



(III.17) 



coinciding with ballistic length b{X). Thus the disordered 
M-stack in the long wavelength region presents a unique 
example of a one-dimensional disordered system in which 
the localization length defined as transmission decrement 
of sufficiently long stack, differs from the reciprocal of the 
Lyapunov exponent. 

The qualitative picture of transmission and localiza- 
tion properties of the symmetric mixed stack described 
above, remains correct in much more general case where 
statistical properties of the r and / layers are different 
and distributions of the fluctuations and thicknesses 
are not rectangular. The only distinction we expect, 
is that localization and crossover lengths will have 
different wavelength dependence that will result in 
more complicated structure of ballistic region like that 
considered below for H-stack (see Section III C| below) . 



10^ 

Jt 10^ 
10^ 



10 



10" 



10" 




FIG. 2: (Ref. |44], color online) Transmission length It vs A 
for M-stack (thick solid line, direct simulation and calcula- 
tions based on WSA recurrence relations) and H-stack (thick 
dashed line, direct simulation). Asymptotics of the localiza- 
tion length / : the short wavelength asymptotic (thin dotted 
line), and the long wavelength asymptotics — thin solid line 
for the M-stack and a thin dashed line for the H-stack. 

To check the WSA theoretical predictions formulated 
above we provided a series of numerical calculations. 
They were made for the lossless stack with uniform 
distributions of the fluctuations S^^\ S^^\ with widths 
of Qjy = 0.25 and Qd = 0.2, respectively and included 
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(a) direc t simu l ations based on the exact recurrence 
relations (11.16), ( IL17| ); (b) the weak scattering analysis 
for the transmission length. In all cases, unless otherwise 
is mentioned, the ensemble averaging is taken over 
Nr = 10^ realizations. 

Throughout this Subsection we considered only M- 
stacks. Nevertheless, to emphasize the main features of 
the transmission in metamaterials, compare transmission 
spectra for a M-stack of = 10^ layers and a H-stack 
of length N = 10^ plotted in the same Fig. [2| Both 
stacks are sufficiently long: for the shortest of them pa- 
rameter NQl is 62.5 ^ 1. There are two major differ- 
ences between the results for these two types of samples: 
first, in the localized regime {N ^ /t), the transmission 
length of the M-stack exceeds or coincides with that of 
the H-stack; second, in the long wavelength region, the 
plot of the transmission length of the M-stack exhibits a 
pronounced bend, or kink, in the interval A G [10^, 10^], 
while there is no such feature in the H-stack results. 

Fig. [2] demonstrates an excellent agreement of ana- 
lytical and numerical results: the curves obtained by 
direct numerical simulations and by calculations based 
on the weak scattering approximation (WSA) are in- 
distinguishable (solid line). The short and long wave- 
length behavior of the transmission length is also in excel- 
lent agreement with the calculated asymptotics in both 
regimes. The characteristic wavelengths of this mixed 
stack are Ai ~ 148 and A2 ~ 839. Therefore, the re- 
gion A < 148 corresponds to localized regime, whereas 
longer wavelengths, A > 839, correspond to the ballis- 
tic regime. Thus the kink observed within the region 
^1 ^ ^ ^ ^2 describes crossover from the localized to 
the ballistic regime. The long wave asymptotic of the 
ballistic length, as we saw below, coincides with that of 
reciprocal Lyapunov exponent. Therefore the difference 
between localization and ballistic lengths of the M-stack 
simultaneously confirms the difference between localiza- 
tion length and reciprocal Lyapunov exponent in local- 
ized regime. 

More detail numerical calculations of transmission 
length, average reflectance, and characteristic wave- 
lengths of the M-stacks with various sizes also demon- 
strate an excellent agreement between direct simulations 
and WSA based calculations thus completely confirming 
the theory presented above^^. 

Until now, we have dealt only with the transmis- 
sion length /t(A), which was defined through an average 
value. However, additional information can be obtained 
from the transmission length /at (A) for a single realiza- 
tion, 



1 ^ In \Tn\ 
In N ' 

In the localized regime, i.e. for a sufficiently long M- 
stack with N ^ the transmission length for a single 
realization /at (A) is practically non-random and coincides 
with /t(A) and while in the ballistic region it fluctuates. 
The data displayed in Fig. |3] enables one to estimate the 
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FIG. 3: (Ref. [44], color online) Transmission lengths It (solid 
black line) and the transmission length for a single realization 
In (dashed blue line) vs A for a M-stack with Qjy = 0.25, Qd — 
0.2 and = 10^ layers. Each separate point corresponds to 
a particular wavelength with its own realization of a random 
stack. 



difference between the transmission length /t(A) (solid 
line) and the transmission length /at (A) for a single ran- 
domly chosen realization (dashed line), and the scale of 
the corresponding fluctuations. Both curves are smooth, 
coincide in the localized region, and differ noticeably in 
the ballistic regime. The separate discrete points in Fig. 
[3] present the values of the transmission length /at (A) cal- 
culated for different randomly chosen realizations. It is 
evident that fluctuations in the ballistic region become 
more pronounced with increasing wavelength. 



C. Homogeneous Stack 

For an H-stack composed entirely of either normal ma- 
terial or metamaterial layers, the transmission length ob- 
tained within the WSA is 



Re 



(a 



1 



1 1 



1 — exp 



H-f)]} 



(IH.18) 



where / is crossover length (HI. 10) and is ballistic 



crossover length defined by equation 



k 



H-stack localization length / is 



-(ln|t|) - Re 



1 - (t^) 



(IH.19) 



where r, t are transmission and reflection coefficients of a 
R layer (for L layers, they should be replaced by r* and 
however this does not change the final result due to 
real part operation Re). 
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Here we consider the simplest lossless model (a = 0) 
with only refractive index disorder (i.e., Qd = 0). In con- 
trast to M-stack case(see Section IIIB below), where a 



minimal model manifesting all common features of the 
M-stack transmission properties necessarily includes ad- 
ditional random parameter (in previous Subsection it is 
layer thickness), for H-stack it is sufficient to include only 
one such parameter. As earlier, we assume uniform dis- 
tribution of refractive index fluctuations with the width 
2Qy. In this case, the short wave asymptotic of the lo- 
calization length coincides with that of M-stack (III. 12) 



and similarly to the M-stack, transmission through short 
H-stacks with N < Q"^ is always ballistic. So below we 
consider long stacks NQl ^ 1. 

In the long wave re gion A ^ 1, the three characteristic 
lengths entering Eq. (III. 18) asymptotically are 



3A2 



47r 



. (III.20) 



The main contribution to the long wave and short wave 
asymptotic of the localization length is related to the first 
term in Eq. (III. 19). Thus, the localization length of the 



H-stack in these two limits is well described by the single 
scattering approximation. The long wave asymptotic of 
the H-stack localization length differs from that of M- 
stack and coincides with that of its reciprocal Lyapunov 
exponent (III. 17) and ballistic length (HI. 15). 



We calculated also H-stack Ly apunov exponent. It is 
described by the same equation (HI. 16) as that for M- 



stack, thus the reciprocal Lyapunov exponents f or bot h 
types of stacks have the same asymptotic form (HI. 17). 



This coincidence was established analytically in a wider 
spectral region in^^. 

Long H-stacks with N ^ Q"^ in the long wave re- 
gion A ^ 1 manifest both ballistic and localized behav- 
ior. Transition between these regimes is gove rned b y 
two characteristic wavelengths defined by Eq. (III.ll). 



Similarly to the M-stack case, they are proportional to 
Qy\fN ^ differ only by a numerical multiplier and satisfy 
the inequality \i{N) < \2{N). 

At starting part of long wave region 1 <C A <C 
Ai(A^) transmission length It coincides with the local- 
ization length / and has asymptotic described by Eq. 
(III.20). Then after passing transition region \i{N) <C 
A <C \2{N) the ballistic regime A ^ \2{N) starts. In 
this regime, transmission length coincides with ballistic 
length h{\) described by equation 



3A2 



6(A) 



NQl 



( ■ UN) \ 

/ Sin \ 



12 



2A 



A3(^') 



2A 



A3 (TV) =47rAr, 



(III.21) 



obtained by expansion of the exponent exp (—TV//) in Eq. 

(IhlTs). 



Due to appearance of additional characteristic wave- 
length Xs{N) determined by equation N = li){Xs{N)) 



where is the ballistic crossover length (III.20), ballistic 



region is naturally divided onto two subregions. The first 
of them defined by inequalities X2{N) <C A <C A3(A^) is 
near ballistic region where ballistic length coincides with 
localization length 



3A2 



(IH.22) 



Thus crossover from localized regime to ballistic one is 
not accompanied by any change of transmission length. 
In the ballist ic tran sition region A ~ A3 (A/"), the second 
term in Eq. (HI. 21) becomes essential leading to oscilla- 



tions of ballistic length. Finally in the far ballistic region 
A Xs{N), expansion of the sine in Eq. (HI. 21) shows 
that for long stacks the second term in this equation dom- 
inates and far ballistic length is 



3A2 



12 18A2 



(III.23) 



The region A > A3 (A'') possesses a simple physical in- 
terpretation. Indeed, in this subregion, the wavelength 
essentially exceeds the stack size and so we may consider 
the stack as a single weakly scattering uniform layer with 
an effective dielectric permittivitj*'*^' 



£eff 



Substitution this value to the text-book formula for re- 
flectivity of the uniform sample lea ds imm ediately to the 
far long wavelength ballistic length (III.23). We note that 



because of the effective uniformity of the H-stack in the 
far ballistic region, the transmission length on a single 
realization is a less fluctuating quantity than that in the 
near ballistic region, where it fluctuate strongly as over 
entire ballistic region for M- stacks. 

Numerical calculations for H-stack, show an excellent 
agreement between direct simulations and calculations 
based on WSA: corresponding curves can not be distin- 
guished. Figure [2] explicitly demonstrates that transmis- 
sion length preserves the same analytical form in local- 
ized long wave region and near ballistic region. For the 
considered stacks with A^ = 10^, the transmission spec- 
trum features corresponding to transition between two 
ballistic subregions can not be manifested. Indeed, the 
transition occurs at the wavelength A ~ 10^ that is out 
of range in this figure. 

To study the crossover from near to far ballistic behav- 
ior, consider the transmission lengths of H-stacks with 
N = 10^ and 10^ over the wavelength range extended 
up to A ~ 10^ plotted in Fig. [ij The transition from 
the localized to the near ballistic regime occurs without 
any change in the analytical dependence of transmission 
length, however the crossover from the near to the far 
ballistic regime is accompanied by a change in the ana- 
lytical dependence that occurs at A = A2(A'), which for 
these stacks is of the order of 10^ and 10^ respectively. 
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FIG. 4: (Ref. [32, color online) Transmission length It vs A 
for H-stacks of = 10^ (solid line), and lO'^ (dotted line) lay- 
ers (numerical simulation and WS A) . Long wave asymptotics 
for the ballistic length in the near and far ballistic regions are 
plotted in thin solid lines. 
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FIG. 5: (Ref. [44], color online) Transmission lengths It (solid 
black line) and the transmission length for a single realization 
In (dashed blue line) vs A for a H-stack with Qu = 0.25, Qd = 
0.2 and N = 10^ layers. Each separate point corresponds to 
a particular wavelength with its own realization of a random 
stack. 



for a single randomly chosen realization (dashed line) are 
plotted. Like the M-stack case, the H-stack single real- 
ization transmission length in the near ballistic region is 
a complicated and irregular function, similar to the well 
known "magneto fingerprints" of magneto-conductance 
of a disordered sample in the weak localization regim^^. 
This statement is supported by displayed in Fig. [5] the 
set of separate discrete points, each of them presenting 
transmission length calculated for a different randomly 
chosen realization. 



The crossover is accompanied by prominent oscillations 
described by Eq. (111.21). Finally, we note that the ver- 



tical displacement between the moderately long and ex- 
tremely long wavelength ballistic asymptotes does not 
depend on wavelength, but grows with the size of the 
stack, according to the law 



In 



In 



NQl 
12 ' 



which stems from Eqs. (III.22), (III.23). 



Detailed study of the average reflectivity of the fi- 
st acks with various lengths at all long wave regioiP^ 
also completely confirm theoretical predictions formu- 
lated above. 

Consider now statistical properties of the H-stack 
transmission length on a given realization /at (A). For very 
long stacks N ^ oo this length becomes practically non- 
random in both localized region due to self- averaging of 
Lyapunov exponent, and far ballistic region because due 
to self- averaging nature of the effective dielectric permit- 
tivity. For less long stacks, transmission length It fluc- 
tuates also even in the far ballistic region, however for 
sufficiently long stacks these fluctuations are essentially 
suppressed since they must vanish in the limit as ^ oo. 
This is demonstrated by Fig. [5] where the transmission 
length It (solid line) and the transmission length /at (A) 



D. Transmission Resonances 

An important signature of the localization regime is 
the presence of transmission resonances (see, for exam- 
ple, Refs.^^^^*^, which appears in sufficiently long, open 
systems and which are the "fingerprints" of a given re- 
alization of disorder. These resonances manifest them- 
selves as narrow peaks of transmittivity \Tn\'^ on a gi ven 
realization as a function of wavelength A. Figure [6p^ 
presents a single realization of the transmittance |T/vp 
as a function of A for a M-stack (dashed line) and for 
the corresponding H-stack of = 10^ layers (solid line). 
It is evident that the resonance properties exhibited by 
homogeneous and mixed media, serve as another (in ad- 
dition to the behavior of the localization length) discrim- 
inating characteristic of these two media. Indeed there 
are no resonances for the M-stack for A > 4, while the 
disordered homogeneous stack exhibits well pronounced 
resonances over the entire spectrum. 

Note that the dotted curve in Fig. [6] describes reso- 
nance properties of periodic = comparatively short 
M-stack with only refractive index disorder (RID). Im- 
portant feature of such a stack is the lack of phase ac- 
cumulation over its total length: in the particular real- 
ization of Fig. [6) the accumulated phase of the wave in 
the mixed stack never exceeds 7r/2. Therefore to subdue 
such a suppression of the phase accumulation, one need 
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FIG. 8: (Ref. |44], color online) Single realization transmit- 
tance |T|^ vs A for M-stack of iV = 10^ layers with Qu = 0.25. 
Solid line corresponds to an M-stack with Qd = 0.2, and 
the dashed line to M-stack with no thickness disorder, i.e., 
Qd = 0.0. 



FIG. 6: (Ref. [43j, color onhne) Transmittance |T|^ vs A for 
a sine: le realization (Q = 0.25, N = 10^). Solid: normal H- 
stack, dotted: M-stack. 



or to enlarge essentially the stack size, or to switch on 
an additional (thickness Qd or magnetic permittivity /i) 
disorder. 
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FIG. 7: (Ref. [44], color online) Single realization transmit- 
tance |T|^ vs wavelength A for RID M-stacks with Qu = 0.25 
and = for AT = 10^ layers (sohd line) and N = 10^ layers 
(dotted line). 

The first possibility is demonstrated in Fig [7| where 
transmittance spectra |Tp(A) for a realization, of two 
different M-stack with two lengths = 10^ and N = 
10^ and only refractive index disorder is displayed. It 
is readily seen that while the resonances in the shorter 
stack (dashed line) at A > 5 do not exist at all, they do 
appear in the same region for the longer sample (solid 
line). 

The second way to generate transmission resonances is 
to introduce additional disorder. This is confirmed by the 



transmittance spectra for a realization, of two M-stacks of 
the same size N = 10^ with only refractive index disorder 
(dashed line), and both (thickness and refractive index) 
types of disorder (solid line), plotted in Fig. [s] It is clear 
that while the RID M-stack with this length, is too short 
to exhibit transmission resonances at A > 3, resonances 
do emerge at longer wavelengths for the M-stack with 
thickness disorder. 

Transmission resonances are responsible for the differ- 
ence between two quantities that characterize the trans- 
mission, namely transmittance logarithm (ln|Tp) and 
logarithm of average transmittance ln(|Tp). The for- 
mer reflects the properties of a typical realization, while 
the latter value is often very sensitive to the existence 
of almost transparent realizations associated with the 
transmission resonances. Moreover, in some cases namely 
small number of such realizations contribute mainly to 
the average transmittance. 

Thus the ratio of the two quantities mentioned above 



(Inm^ 
ln(|TP) • 

is a natural characteristic of the transmission resonances. 
In the absence of resonances, this value is close to unity, 
while in the localization regime s > 1. In particular, 
this ratio takes the value 4 in the high-energy part of the 
spectrum of a disordered system with Gaussian white- 
noise potentiaP. 

Consider the ratio 5(A) as a function of the wavelength 
for RID M- and H-stacks and for the corresponding stacks 
with thickness disorder plotted in Fig. [9) In all cases, the 
stack length is N = 10^. It is evident that for the RID 
M-stack s{X) ~ 1, i.e.. the length of this M-stack is too 
short for the localization regime to be realized. In other 
three cases, however, s(A) > 2, which means that the 
localization takes place even in such comparatively short 
stack. 
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FIG. 9: (Ref. 44 , color online) Ratio s(A) vs. wavelength A 
for Q^J — 0.25 and the stack length N — 10^. Solid and dashed 
curves are for the RID H- stack and H-stack with Qd = 0.2, 
respectively. The middle dashed-dotted curve is for an M- 
stack with Qa = 0.25, and the bottom dotted line is for a 
RID M-stack. 



Characteristic angle conserves its direct geometrical 
meaning for incidence angle < Oc { subcritical incidence 
angle) where critical angle is 



For the supercritical incidence angle 9 > Oc^ the values 
of are complex. 

Below we mention only final asymptotical expressions 
for some characteristic lengths of the problem in the 
typical cases. We take into account both two types of 
disorder however in all final results we keep only the 
leading terms and omit the higher order corrections with 
respect to the refractive index and thickness fluctuations 
Qu,d' 

In the short wave limit, localization length is the same 
for M- and H-stacks. In the subcritical region of incidence 
angles it is 



E. Polarization Effects 

The results obtained above for normal incidence can 
be easily generalized to the case of oblique incidence. 
Here all characteristic lengths and wavelengths depend 
on angle of incidence and s- and p-polarizations should 
be considered separately. Qualitatively new features ap- 
pear: essential enlightening in the vicinity of Brewster 
angle and appearance of supercritical regime induced by 
total internal reflection. We describe these new proper- 
ties within the frameworks of the model defined in the 
previous Subsection |III A 



General expressions for transmission length for both 
M-stac ks (Eqs. ( |III.7D - ( |III.10D ) and H-stacks (Eqs 
(|III.18|) - ( III.19[ )) as well as general expressions (III. 2) 



for transmission and reflection coefficients of a single 
layer remain the same as in the case of normal incidence. 
However explicit expressions for the parameters entering 
these coefficients are changed. Fresnel interface refiection 
coefficient is now given by 



(III.24) 



Z cos 6iy — cos 
^ Z cos + cos d ' 
5 — polarization 

Z, p — polarization. 

Here characteristic angle and the layer impedance Z 
relative to the background (free space) according Eq. 
UlLil are 



sm 



cos6>i, = y 
Then the phase shift /3 is now 



\i 1 

7 ~ iT^' 



P = kducosOjy^ k = 27r/X. 



(in.25) 



Ql 



12cos4|9 



1 s — polarization, 
cos^ 26 p — polarization. 



Note that for p-polarization, this expression acquires an- 
gle dependent multiplier that vanishes at the Brewster 
angle = 7r/4. Accounting for the next term we obtain 
the localization length at the Brewster angle 



1 = ^ 

which is times larger than that far from the Brew- 
ster angle and than that for 5— polarization in the same 
shortwave limit. 

At the incidence angle > Oc, total internal refiection 
occurs and the WSA fails. If the supercriticality 6 — Oc is 
not extremely small, then the exponent 2iP in Eq. (HI. 2) 



is real and negative and thus the magnitude of the single 
layer transmission coefficient is exponentially small. This 
results in the attenuation length for both polarizations 



1 

^att 

k sin^ 



Im(/3) = k lm{dV sin^ — v'^) 



8Q. 



(tt — 2^0 — sin26>o), sin6>o = 



smt 



smf 



Due to (X /c dependence, in the short wave limit latt oo 
and transmission length in supercritical region of the an- 
gles of incidence coincides with the attenuation length. 
However for the same reason at long waves attenuation 
contribution can be neglected and the main contribution 
to the transmission length is due to Anderson localiza- 
tion. 

In the long wave region, H- and M-stacks demonstrate 
different behavior and we describe them separately. 
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a) Homogeneous stacks. 

For 5-polarization, long wave asymptotic of the transmis- 
sion length is similar to that for normal incidence (III. 21 ) 



3A2 cos2 



NQl 
12 



2ttN cos e ' 



sm ■ 



A 



V 



271 N cos 



X 



This expression describes localized region as well as both 
ballistic subregions. 

In the case of p— polarized wave, the localization length 
is given by 



1 2n^Qlcos^20 ^2g4 
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At Brewster angle - 
transmission length is 



the first term vanishes and 



1 ^ Wn^Qt 
I ~ 45A2 • 



(III.26) 



b) Mixed stacks. 

Reciprocal transmission length for s-polarized wave is 
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FIG. 10: (Ref. [45 , color online) Transmission length It ver- 
sus A for a M- stack in p-polarized light with Qi, = 0.1, 
Qd = 0.2 and AT = 10^, at the Brewster angle = 45° red solid 
line. The blue dashed line shows results for 5-polarization and 
a H-stack, re-plotted for comparison. 



k^Qi cos^ 29 

3 cos^ 9 
1 - f{Na^) 



1 

2 ~ 2 + cos2 26' + Ccos4 6' 

-^^(2 + cos2 20 + Ccos4 
3 cos^ 



, (III.28) 



The localization length is deduced from Eq. (III. 28) by 
taking the limit as ^ oo 



1 



3 cos^ 
3 cos^ 6 



3 + C cos4 (9 
(3 + Ccos'^ 0). 



(III.27) 



where the function / and parameter ( are defi ned in Eqs. 
(III. 8) and (IIIB) correspondingly. Equation (III. 27) de- 



scribes the transition from localization to ballistic prop- 
agation at long wavelengths. In the limit N ^ oo trans- 
mission length tends to localization length 



3A^ cos^ 6 3 + Ccos^<9 
~ 27r2Q2 i + (cos4(9' 

while the opposite extreme, i.e., as TV ^ 0, gives the 
ballistic length 



3A^ cos^ 

which coincides with that for a H-stack in 5-polarization. 

For p-polarized waves incident at angles away from the 
Brewster angle, the transmission length is given by: 



/ = 



3A2 cos^ (9 2 + cos^ 2(9 + C cos"^ (9 
27r2 Ql cos2 2(9 cos^ 2(9 + C cos^ ' 



Correspondingly, the ballistic length is obtained by cal- 
culating the limit as ^ 



b = 



3A^ cos^ 
27r2Q2cos2 2(9' 



At the Brewster angle = 7r/4, accounting for the 
higher order corrections to r.h.s. of Eq. ( III. 28 ) we obtain 
the transmission length the same result (III. 26) that for 
H-stack. 

All analytical predictions are confirmed by numerical 
calculations. As in the case of normal incidence the- 
oretical curves based on WSA approximations mostly 
can not be distinguished from those obtained by direct 
simulations. The results obtained mostly similar to 
those of normal incidence. Therefore here we mention 
only some of them which differ from presented above. 



In Fig. 



10 the transmission length spectrum of an 
M-stack of length N — 10^ in p-polarized light with 
other parameters Qy — 0.1, Qd = 0.2, Nr = 10^, 
and the incidence angle = 7r/4 is displayed. The 
chosen angle of incidence is less than the critical angle 
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FIG. 11: (Ref. |35], color online) Transmission length It 
versus A for a M-stack in s-polarized light with Qi, — 0.1, 
Qd — 0.2 and N — 10^, and for the supercritical incidence 
angle = 75°. Red solid curve: numerical simulations; Blue 
dash curve: analytic form. 



Oc = arcsinO.9 = 64.16° and coincides with the Brewster 
angle for the single layer with mean refractive index 
= ±1. The results of the numerical simulation and 
the WSA analytical forms coincide and are displayed 
by a single red solid line. Localization occurs for 
A < Ai ~ 19, while the transition from localization to 
ballistic propagation occurs at A ~ Ai. In contrast to the 
case of 5-polarization, this transition is not accompanied 
by a change of scale and is given by the same wavelength 
dependence. Transition from near to far ballistic length 
is accompanied by oscillations of transmission length 
which are much more pronounced in comparison to the 
case of normal incidence. 

Consider now a supercritical case where the angle of 
incidence 6 = 75° exceeds the critical angle. In Fig. 
pT] we present the transmission length spectrum for 5- 
polarized light is presented. The results of both the ex- 
act numerical calculation (red solid line) and the ana- 
lytic form (long dashed blue curve) are displayed. The 
short wave (dashed dotted line) and the long wave (black 
dashed line) asymptotic of the transmission length, re- 
spectively coincide with the numerical results for A < 1 
and 200 < A. In the intermediate region 1 < A < 200, 
however, the theoretical description underestimates the 
actual transmission length since the WSA is no longer 
valid for the chosen, supercritical angle of incidence. For 
p-polarization, the results are qualitatively the same, but 
with the discrepancy at the intermediate wavelengths 
even more pronounced. 

We consider also the angular dependence of the trans- 
mission length for mixed stacks. In Fig. 12 the trans- 



mission length It as a function of the angle 6 for a stack 
of length N = 10^ at the two wavelengths A = 0.1 and 
A = 1 is displayed. In either case, the calculated trans- 
mission length does not exceed the stack length and so, 
for subcritical angles, our calculations display the true 
localization length. For the shorter wavelength A = 0.1, 
the form of the transmission length for both polarizations 



It 10^ 




FIG. 12: (Ref. |45^, color online) Transmission length It 
versus incidence angle for a mixed stack with Qi, = 0.1, 
= 0.2, for (a) A = 0.1 (upper panel), and (b) A = 1 (lower 
panel) . The top and bottom curves are respectively for p- and 
5-polarizations. 



is similar to that observed for homogeneous stacks. 
Fig. 12 ^b) displays results for an intermediate wave- 
1 with the lower solid red and blue dashed 



length A 

curves respectively displaying the results of numerical 
simulations and analytical predictions for 5-polarization, 
(bottom curves), while the upper solid green and brown 
dashed curves display simulations and analytical predic- 
tions for p-polarization. The agreement between simula- 
tions and the theoretical form is again excellent for an- 
gles of incidence less then the critical angle, 6 < 6c, while 
for angles greater then the critical angle, the discrepan- 
cies that are evident are again explicable by the breaking 
down of the WSA at extreme angles of incidence. 



F. Dispersive Metamaterials 

Real metamaterials always are dispersive materials. 
Here we consider a dispersive model of the stack com- 
posed of metalayers with the same thickness d and ran- 
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dom dielectric permittivity and the magnetic permeabil- 
ity described by Lorentz oscillator model 



£(/) 



f2 _ f2 

J ep J e 

f2 _ f2 

J mp J m 



(111.29) 
(111.30) 



Here / is circular frequency, fm and fe are the resonance 
frequencies and 7 is the phenomenological absorption 
parameter. In this model, disorder enters the problem 
through random resonance frequencies so that 



fe /e(l + ^e), fm = /m(l + ^m), 

where /e,m = {fe,m) are the mean resonance frequencies 
(with the angle brackets denoting ensemble averaging) 
and ^e,m are independent random values distributed uni- 
formly in the ranges [— Qe,m, Qe,m]- The characteristic 
frequencies f^p and f^p are non-random. Therefore, in 
lossless media (7 = 0), both the magnetic permeability 
and the dielectric permittivity vanish with their mean 
values, e{f) = {s{f)) and fl{f) = (/i(/)), at frequencies 
/ = fep and / = fmp respectively, i.e.. 



f^ifmp) = l^ifmp) 0, e{fep) 



I ep) 



0. 



Following Ref.^^^^, in our numerical calculations we 
choose the layer thickness d = 0.003m and the values of 
characteristic frequencies f^np = 10. 95 GHz, fmQ = fm = 
10.05GHz, fep = 12.8GHz, feo = fe = 10.3GHz, and 
7 = lOMHz, which fit the experimental data given in 
Ref.!^. That is, we are using a model based on experi- 
mentally measured values for the metamaterial parame- 
ters. Then we choose the maximal widths of the distribu- 
tions of the random parameters ^e,m as Qe.m — 

5 X lO-'^ 

corresponding to weak disorder. 

We focus our study on the frequency region 
10.40GHz < / < ll.OOGHz. In the absence of absorp- 
tion and disorder, for these frequencies the dielectric per- 
mittivity and the magnetic permeability of the metama- 
terial layers vary over the intervals —26.9 < e < —2.9 
and —1.64 < ja < 0.055. The refractive index is neg- 
ative in the frequency range 10.40GHz < / < fmp = 
10. 95 GHz, as shown in the inset of Figjl3j However, at 
fmp = 10.95GHz, the magnetic permeability changes sign 
and the metamaterial changes from being double nega- 
tive (DNM) to single negative (SNM). As we show later, 
such changes have a profound effect on the localization 
properties. 

We study the transmission of a plane wave either s - 
or p - polarized and incident on a random stack from free 
space with an angle of incidence ^o- 

In the previous Subsections, we have described and 
used an effective WSA method developed and elaborated 
in Refs^^ for studying the transport and localization 
in random stacks composed of the weakly reflecting lay- 
ers. 



In the dispersive case, the reflection from a single layer 
located in free space is not necessarily weak, in which 
instance the method seems inapplicable. However, we 
can replace each layer by the same layer surrounded by 
infinitesimally thin layers of a background medium with 
permittivity and permeability given by the mean values 
of e{f) = {s{f)) and /i(/) = (m(/)) respectively. In the 
considered case of weakly disordered stacks, we can use 
the WSA approximation for all layers beside two " leads" 
connecting the stack with free space from the very left 
and the very right its ends. The localization characteris- 
tics which are intrinsic properties of the stack do not feel 
the leads. Their role is restricted by only change the cou- 
pling conditions to the random stack through the angle 
of incidence transforming it from its given value ^0 out- 
side the lead to the frequency dependent refracted value 
6>5 inside the lea d. These angles are related by Snell law 
sin6>o = sin \/s{f)jl{f). It is important to note that, 
while in the localized regime the input and output leads 
are of no significance, they do play a crucial role when 
localization breaks down (see below). 



The single layer scattering is described by Eqs. (HI. 2) 
where according Eqs. (III.l) and(III.25) 



Pn kdVn COS On, T^n = V^n/^n, 



(IH.31) 



27r 27r/ 

and /c = — - = is the free space wave number. The 

A c 

interface Fresnel reflection coefficient pn is given by 



Pn 



Z5 COS Ob — Zn COS On 



Z\) cos 6>5 + Zn COS On ' 

The impedances Z^, and Zn are 

pie p-polarization, 
sIp s-polarization, 

\/ Pn/^n p-polarization. 



(IH.32) 



\/ Enl Pn s-polarization. ' 



and the angles 6^5 and On satisfy Snell's law 



V sm 6/5 



sm e/o. 



£^.(IH.33) 



General WSA expressions (HI. 19) and (HI. 9) for lo- 



calization length of mono-type and mixed stacks remain 
valid for the stacks composed of dispersive stacks. To 
study localization properties of such stacks we should 
insert there the same single layer scattering coefficients 
(HI. 2) wi th disp ersive phase shift (HI. 31) and Fresnel co- 
efficient ( [IIL32| . 



Dispersion affects essentially the transport proper- 
ties of the disordered medium. In particular, it can 
lead to suppression of the localization either at some 
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angle of incidence, or at a selected frequency, or even 
in a finite frequency range. Below we consider the 
two first cases for the H-stack composed of L-layers. 
The third case will be considered further in Section llll Gl 



ikL 1 ^ , 



In the presence of dispersion, the long-wave asymptotic 
of the localization length is 



In a sufficiently long stack, £ ^ ^ikLe and the transmit- 
tance T= |7ii|~^ is given by 



1_ TT^d^ , {e')-{e)' 

I 2X^f)\ + (e)2 



(III.34) 



where n and s are given by Eqs. (III. 30), (111.29), and 
frequency-dependent wavelength in the medium 



A(/) 



fV^fUf) 



and can be large even when the wavelength of the incident 

c 

signal, A = — , is small. 

Accordingly, the inverse localization length 



r 



fe{f)^{f) 



becomes small not only at low frequencies / ^ but 
also in the vicinity of ja- or e- zero points. For example, 
as the frequency approaches the /i-zero point from be- 
low, i.e., / f~p, in a H-stack of metamaterial layers, 
for any realization, is proportional to the differ- 
ence {fmp — f) and the expression for localization length 
diverges as {fmp — /)~^- Formally, this divergence can 
be treated as delocalization, however the limiting value 
1// = means nothing but the absence of exponential 
localization. Moreover, when the localization length be- 
comes larger than the size of the stack, ballistic trans- 
port occurs and the transmission coefficient is determined 
by transmission length, rather than by the localization 
length. 

To calculate the transmission coefficient for this case 
we consider, for the sake of simplicity, a stack with only 
e-disorder. Here the transfer matrix of the n-th layer at 
/ = fmp has the form 



Tn T{en) 



1 + en en 
-en 1 - er: 



where = ikdsn/'^. 

As a consequence of the easily verified group property 



r(ei)r(62) =r(ei+e2), 
it follows that the stack transfer matrix T is 



N 



1 + f £ 
-£ l-£ 



fkLe{f)\ 



1 



Thus, at the frequency fmp, the transmittance of the 
sample is not an exponentially decreasing function of 
the length L (as is typical for ID Anderson localiza- 
tion). It decreases much more slowly, namely, accord- 
ing to the power law T ex The explanation of such 
a decrease is that at a /i-zero point (/ = fmp), the re- 
fractive index vanishes together with the phase shift 
Pn = kdun COS On across the layer, thereby destroying the 
interference, which is the main cause of localization. An- 
other form of the explanation is that the effective wave- 
length inside the stack tends to infinity when /i ^ and 
exceeds the stack length. Obviously, such a wave is in- 
sensitive to disorder and therefore cannot be localized. 

In the limit as the frequency approaches the /i-zero 
frequency, from above, i.e., / Z^^, the medium is 
single-negative and e/i < 0. For frequencies / not too 
close to fmp, the radiation decays exponentially inside the 
sample due to tunneling, and in the absence of dissipation 
the decay rate is: 



1 



kd^-{li){e) 



(IH.35) 



where 



Thus, as we approach the /i-zero frequency from the 
right, the formally-calculated localization length diverges 
as / (X (/ — fmp)~^^'^ i-e. much more slowly than for the 
left-hand limit for which / cx {fmp — /)~^- The trans- 
port properties in the vicinity of the e-zero frequency fep 
can be considered in a similar manner. Waves are also 
delocalized in the more exotic case when both dielectric 
permittivity and magnetic permeability vanish simulta- 
neously. The vanishing of both ji and e simultaneously 
can happen at Dirac points in photonic crystalJ^. 

The use of off-axis incidence from free space for fre- 
quencies for which ji oi e are zero is not an appropriate 
mechanism for probing the suppression of localization. 
In such circumstances, tunneling occurs and the localiza- 
tion properties of the stack are not "accessible" from free 
space. Nevertheless, suppression of localization can be 
revealed using an internal probe, e.g., by placing a plane 
wave source inside the stack, or by studying the corre- 
sponding Lyapunov exponent. Both approaches show to- 
tal suppression of localization at the frequencies at which 
dielectric permittivity or magnetic permeability vanish. 

In such circumstances, each layer which is embedded 
in a homogeneous medium with material constants given 
by the average values of the dielectric permittivity and 
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magnetic permeability, is completely transparent, with 
this manifesting the complete suppression of localization. 
However the "delocalized" states at the zero-/i or zeio-e 
frequencies are in a sense trivial, corresponding to fields 
which do not change along the direction normal to the 
layers. 

Another example of the suppression of localization is 
related to the Brewster anomaly. As we saw above, in 
a non-dispersive mixed stack with only thickness disor- 
der, delocalization of p-polarized radiation occurs at the 
Brewster a ngle of incidence. At this angle, the Fresnel co- 
efficient p (III. 24) and, therefore, the reflection coeflicient 



(III. 2 ) as well, vanish for any frequency, thus making each 



layer completely transparent. 

In the presence of dispersion, the same condition p = 
leads to more intriguing results. In this instance, 
frequency-dependent angles, at which a layer becomes 
transparent, exist not only for p-polarization, but also 
for an 5-polarized wave. This means that the Brewster 
anomaly occurs for both polarizations, with the corre- 
sponding angles. Op and ^g, being determined by the con- 
ditions 



Note that this equation may be satisfied at multiple fre- 
quencies depending on the form of the dispersion. 

The case of only magnetic permeability disorder, e = 
is described by similar equations which are obtained by 
replacement s ^ p in Eqs. (??) - (III.40). 



For disorder in both the permeability and the permit- 
tivity, the existence of a Brewster anomaly a ngle de - 
pends, in accordance with Eqs. ( III.36| ) and (III.37), 
on the sign of the quantity ^ = {ep — ep)/{ep — ep). 
If ^ > 0, the Brewster angle exists for 5-polarization, 
while if ^ < 0, it exists for p-polarization. In the case 
^ = 0, the layer and the medium in which it is embedded 
are impedance matched, and thus the layer is completely 
transparent. 

The features of transmission length mentioned above 
are completely confirmed by numerical calculations. 
Consider first the case of normal incidence on a stack 
of = 10^ layers, in which we randomize only the di- 
electric permittivity {Qm = 0) with Qe = 0.5 x 10~^. 
In Fig. [13] the transmission length It as a function of 
frequency / is displayed. The upper curves present the 
lossless case, while the lower curves show the effects of 
absorption (see^^ for details). 



tan^ On 



tan^ Oo — 



e{ep — ep) 
e{ep — ep)^ 
p{ep — ep) 
p{ep — pe) 



(III.36) 
(III.37) 



The right hand sides of these equations always h ave op- 
posite signs. Therefore from Brewster conditions (III. 36) 



and (III. 37) one can find either the Brewster angle and 



corresponding polarization for a given frequency, or the 
Brewster frequency and corresponding polarization for a 
given angle of incidence. 

While, for a stack with only thickness disorder, the 
condition p = can be satisfied for all layers simultane- 
ous ly, whe n e and/or p fiuctuate, the conditions (III. 36) 
or (III. 37) define the frequency-dependent Brewster an- 



gles which are slightly different for different layers. These 
angles occupy an interval within which the stack is not 
completely transparent, but has anomalously large trans- 
mission lengths^^ . 

When only the dielectric permittivity is disordered and 
/i = /i, the Brewster conditions (III.36), (III. 37) simplify 
to 



tan^ Og = 
tan^ Onj = 



1. 



(III.38) 
(III.39) 



In this case, the Brewster condition is satisfied only for 
p-polarization. For weak disorder, the Brewster angle of 
incidence from the effective medium is Op ~ 7r/4. For a 
given frequency /, angle of incide nce fro m free space, ^o, 
should be found from Snell's law (III.33), and for a given 
6>o, the Brewster frequency fp follows from 



kfpMfp) = ^ = V2smdo. (III.40) 

sm Un 
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FIG. 13: (Ref. [46l, color online) Transmission length It vs 
frequency / at normal incidence (Oa = 0) for a metamaterial 
stack without absorption (top curve) and in the presence of 
the absorption (bottom curves). Red solid curves display nu- 
merical simulations while blue dashed curves show the analyt- 
ical predictions. Inset: the real (red solid line) and imaginary 
(green dashed line) part of the metamaterial layer refractive 
index. 

The red, solid curves and the blue, dashed curves dis- 
play results from numerical simulations and the WSA 
theoretical prediction respectively. The top curves rep- 
resent the genuine localization length for all frequencies 



except those in the vicinity of f ^ fn 



10.95GHz 



where the transmission length dramatically increases. 

In the absence of absorption, for frequencies / > 
10.95GHz, the metamaterial transforms from being dou- 
ble negative to single negative (see inset in Fig. 13). 



The refractive index of the metamaterial layer changes 
from being real to being pure imaginary, the random 
stack becomes opaque, and the transmission length sub- 
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stantially decreases. Such a drastic change in the trans- 
mission length (by a factor of 10^) might be able to ex- 
ploited in a frequency controlled optical switch. Across 
the frequency interval 10 AGHz < f < ll.OGHz, theo- 
retical results are in an excellent agreement with those 
of direct simulation. Moreover, for all frequencies ex- 
cept in the region lOAGhz < f < 10.5GHz, the sin- 
gle scattering approximation excellently describes the 
It be havior. Qui te surp risingly, the asymptotic equa- 
tions ( |III.34 ) and (III. 35) are in the excellent agreement 
with the numerical results even over the frequency range 
10.9Ghz < / < ll.OGhz, including in the near vicinity of 
the frequency fmp = 10.95GHz at which /i vanishes. 

Absorption substantially influences the transmission 
length (the lower curve in Fig 13)^^ and smoothes the 
non-monotonic behavior of the transmission length for 
/ < 10.5GHz. The smah dip at / ^ 10.45GHz correlates 
with the corresponding dip in the transmission length 
in the absence of absorption. The most prominent ef- 
fect of absorption occurs for frequencies just below the 
/i-zero frequency fmp = 10.95GHz. While in the absence 
of absorption, the stack is nearly transparent in this re- 
gion, turning on the absorption reduces the transmission 
length by a factor of 10^-10^ for / > 10.7GHz. In con- 
trast, for frequencies / > 10.95GHz, the transmission 
lengths in the presence and absence of absorption are 
nearly identical because here the stack is already opaque 
and its transmittance is not much affected by an addi- 
tional small amount of absorption. 

The case where both disorders of the dielectric permit- 
tivity and magnetic permeability are present, is qualita- 
tively similar to that of the single disorder case considered 
above. 
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FIG. 14: (Ref. [46], color online) Transmission length It vs 
frequency / for Oa = 30° for a metamaterial stack: without 
absorption, p-polarization (top curves), s polarization (mid- 
dle curves); in the presence of absorption (bottom curves). 

In the case of oblique incidence, polarization effects 
become important. In FigfM} the transmission length 
frequency spectrum is displayed for the same metama- 
terial H-stack with only dielectric permittivity disorder 
for the angle of incidence = 30°. Here for frequencies 
/ < 10.55GHz, the transmission length is largely inde- 



pendent of the polarization. Moreover it does not differ 
from that for normal incidence (compare with the top 
curve in Fig 13 ). This is due to the high values of the 
refractive indices at these frequenci es (li^nl > 4), result- 



ing in almost zero refraction angles (HI. 33) for angles of 
incidence that are not too large. 

The transmission length manifests a sharp maximum 
at an angle close to the Brewster angle, as commented 
upon in Ref^^^EH. This is indeed apparent in Fig. 14 for 
the frequency / ^ 10.85GHz. Because only e fluctuates, 
the Bre wster condition is satisfied only for p-polarization 
[III. 38 ) at a single frequency fp ~ 10. 85 2 GHz. The intro- 
uction of additional permeability disorder (not shown) 
reduces the maximum value of the localization length by 
two orders of magnitude. 



Comparison of Figs [13] and [14] shows that the frequency 
of the maximal suppression of localization decreases as 
the angle of incidence increases. At normal incidence it 
coincides with the /i-zero frequency fmp while for oblique 
incidence at Oq = 30° it coincides with the Brewster fre- 
quency fp for p-polarization. 

Absorption strongly diminishes the transmission pro- 
viding the main contribution to the transmission length 
while the permittivity disorder has little influence on the 
transmission length. In this case, the results for both two 
polarizations are therefore practically indistinguishable. 

The transmission properties of a stack with only mag- 
netic permeability disorder at oblique incidence, are sim- 
ilar to those for the case of only dielectric permittivity 
disorder. The key difference is that there is a Brewster 
anomaly for 5-polarization while for p-polarization it is 
absent. 

We consider also the dependence of the transmission 
length on the angle of incidence at a fixed frequency. 
The results for both polarizations are displayed in Fig. 
\TE\ Here we have plotted the transmission length of 
the stack with only dielectric permittivity disorder with 
Qe = 0.5 X 10-2 at the frequency / = 10.90GHz. The 
upper and middle curves in this figure correspond to 
the results for p- and s-polarized waves respectively in 
the lossless case. For s-polarized light, the transmission 
length decreases monotonically with increasing angle of 
incidence, while for p- polarized wave it increases with 
increasing angle of incidence. Such behavior reflects the 
existence of a Brewster angle for p-polarization at the 
Brewster angle = 20°. The red solid curve shows the 
results of simulations, while the blue dashed line is the 
analytic prediction. 

As in the previous cases, in the presence of absorption, 
the results for both polarizations are almost identical (the 
lower curves in Fig. 15). For angles < 30°, the trans- 



mission length is dominated by absorption, while for an- 
gles ^0 > 30° tunneling is the dominant mechanism. The 
results for permeability disorder are very similar to those 
for permittivity disorder. 

For normal H-stacks, the transmission length mani- 
fests exactly the same behavior as for H-stacks comprised 
of metamaterial layers. 
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FIG. 15: (Ref. |36], color online) Transmission length It vs 
angle of incidence for a homogenous metamaterial stack at 
/ = 10.7GHz with permittivity disorder: in the absence of 
absorption — upper curve and for p polarization; middle curve 
is for s polarization and in the presence of absorption and for 
both polarizations, lower curves. 



cancelation of the phase that occurs in alternating L- and 
i^-layers'^^. 




-2 2 

10 10 10 



G. Anomalous Suppression of Localization 

In this Section, we consider the stacks with only re- 
fractive index disorder (RID) i.e. the stacks with = 

= 0. In this limit, there is nothing special for H- stacks. 
Their transmission length demonstrates qualitatively and 
quantitatively the same behavior as was observed in the 
presence of both refractive index and thickness disorder. 
Corresponding formulae for the transmission, localiza- 
tion, and ballistic lengths can be obtained from the gen- 
eral case by taking the limit as Qd 0. 

In the case of M-stacks, however, the situation changes 
markedly. Here suppression of localization in the long 
wave region becomes anomalously large enhancing trans- 
mission length on some orders of magnitude and even 
changing its functional dependence on the wavelengttP^. 
Instead of the universal cx dependence, long wave 
asymptotic of both localization length / and reciprocal 
Lyapunov exponent follows a power law cx with 
much larger exponent m. 

Let us start with some numerical results demonstrat- 
ing such an anomalous growth of the long wave local- 
ization lengths /, of the minimally disordered M-stack 
with only RID. In Fig 16 localization length for M- 



stack with Q = 0.25 is plotted. Solid line in Fig. [16 
corresponds to for the propagation in a M-stack and a 
single realization of = 10^ layers, while the dashed line 
is for the corresponding H-stack with the same parame- 
ters. Within the localization region 1^{X) < 10^, M-stack 
reciprocal Lyapunov exponent grows in the long wave re- 
gion essentially faster than that of H-stack. While for 
H-stack is described by standard exponent m = 2, its 
value for M-stack was estimated as m = 6 and the phe- 
nomenon itself was named as A^ anomaly. The observed 
anomalous suppression of localization was attributed to 
a lack of phase accumulation over the sample, due to the 



FIG. 16: (Ref. |43], color online) Localization length vs. 
wavelength A for Q = 0.25 and N — 10^ layers; solid line is 
for the M-stack, while the dashed line is for the corresponding 
(normal) H-stack. 

Anomalous suppression of localization is manifested 
also in the case of oblique incidence. The next Figure 17 
displays transmission length spectra for a M-stack with 
only refractive index disorder for an angle of incidence 
of ^ = 30°. There is a striking difference between the 
two polarizations: in the case of p-polarized light, there 
is strong localization at long wavelengths (A < 10^), with 
the localization length showing oc A^ dependence. In con- 
trast, the localization length for 5-polarized light is much 
larger and is estimated as ~ A^ dependence as occurs for 
normal incidence. Note that for 5-polarization, anoma- 
lous enlightening manifests itself only in localization re- 
gions in Fig. [TT] which are bounded from above by the 
wavelength limits A < 5,9, and 12 for stacks of length 
TV = 10^ 10^, and 8 x 10^ respectively. 

This asymmetry between the polarizations suggests 
that the suppression of localization is due not only to 
the suppression of the phase accumulation but also to 
the vector nature of the electromagnetic wave. Because 
of the symmetry of Maxwell's equations between the 
electric and magnetic fields, it is to be expected that 
for a model in which there is disorder in the magnetic 
permeability (with £ = ±1) the situation will be inverted 
with anomalous enlightening for p-polarized waves and 
with 5-polarization showing strong localization. 

The results of calculations^^ provided for much longer 
stacks up to = 10^^ qualitatively completely coincided 
with the previous ones. However more detailed studies 
quantitatively occurred slightly different. Generation of a 
least squares fitting It = AX^ to the transmission length 
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FIG. 17: (Ref. 45 , color online) Transmission length It ver- 
sus A for a M-stack with = 0.25, = and 6* = 30° for 
]9-polarized light (cyan dashed dotted curve, N — 10^) and s- 
polarized light (red solid curve, N — 10^; green dashed curve, 
N = 10'^; blue dotted curve, AT = 8 x 10^). 



data, led to a bit surprising conclusions. The best fits 
were m ^ 6.25 for N 10^, m ^ 7.38 for N 10^, 
and even m ^ 8.78, for TV = 10^^. This shows that the 
question about a genuine value of exponent m remains 
still open. 

Consider now the long wave behavior of the localiza- 
tion length in the presence of dispersion. In the panel 
the transmission length spectrum is 



a) of the Fig 18 



plotted in the case of normal incidence, for a small per- 
mittivity disorder of Qe = 0-5 x 10~^. One can im- 
mediately observe significant (up to four orders of mag- 
nitude) suppression of localization in the frequency re- 
gion 10.50GHz < / < 10.68GHz. However, this sup- 
pression seems to have nothing common with observed 
above anomalous enlightening. Indeed, in this case the 
localization length grows with increasing frequency, while 
in the previous studieJ^^^EIl, similar growth has been 
observed with increasing incident wavelength. This is 
demonstrated in FigfT8|3 where the same transmission 
length spectrum is plotted as a function of free space 
wavelength. Thus, the localization length decreases by 
four orders of magnitude, manifesting as an enhance- 
ment, rather than the suppression, of localization with 
increasing wavelength. 

Although at the first sight these findings are in sharp 
contrast with the previous ones, they are correct and 
physically meaningful. In the model studied earlieiE^'EIl, 
the wavelength of the incident radiation largely coincided 
with the wavelength inside each layer. In dispersive 
medium considered here, these two wavelengths differ 
substantially. Accordingly, in FigjlSfc, we plot the 



transmission length as a function of wavelength within 
the stack and obtain results which are very similar to 
those in Refs.l^^M^. To emphasize this similarity, we 
have plotted the transmission length spectrum for three 
different stack lengths: N = 10^,10^,10'^. It is easily 
seen that the suppression of localization in the dispersive 
media is qualitatively and quantitatively similar to 
that predicted in Ref.^^. Corresponding exponent m of 
anomalous enlightening estimated with the help of these 
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FIG. 18: (Ref. [46 , color online) (a) Transmission length It 
vs frequency / for a mixed stack with N = 10^ layers (top 
dotted blue curve), and only dielectric permittivity disorder. 
The bottom curves on all panels (a,b,c) are for a stack with 
N = 10^ layers with both permittivity and permeability dis- 
order (cyan, solid curve display simulation results while the 
dashed, black curve is for the analytic prediction; (b) is the 
same as in (a) but plotted as a function of the free space wave- 
length Ao while on panel (c) we plot the transmission length 
as a function of the averaged wavelength inside of the stack 
normalized to the thickness of the layer, for N = 10^ layers 
(blue dotted top curve), N — 10^ layers (dashed green curve) 
and for N — 10^ layers (red solid curve) respectively. 



results, is m ~ 8.2. 

Enhanced suppression of localization exists in the 
strictly periodic alternative M-stacks with a constant 
layer thickness and only refractive index disorder. By 
other words, in mixed stacks having constant layer thick- 
ness, the dielectric permittivity disorder alone is not suf- 
ficiently strong to localize low-frequency radiation by a 
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standard way. There are many ways to violate these con- 
ditions. It is possible to add thicknesses fluctuations^^ 
or magnetic permeability fluctuation^^, or to introduce 
a weak difference between two constant thicknesses of 
R- and L-layers, or not to change any parameter but 
rearrange randomly the same numbers N/2 of R- and 
L-layers . Each such a violation immediately destroys 
anomalous suppression of localization and restores stan- 
dard long wave asymptotic / ex A^. 

The analytical results obtained above in Section |III[ 
survive in the 5^ ^ limit and predict / oc asymptotic. 
However more detailed investigation shows that WA in 
its form (|II.18|), (|II.19|) fails in this limit 



As was mentioned above, localization length man- 
ifests qualitatively the same behavior as transmission 
length It- At the same time, its calculation is simpler 
than that of It- Lyapunov exponent in minimally dis- 
ordered M-stacks was calculated in^^ using some version 
of the method described in Refs. ^^J and at the 

end of Section |IIB[ The remaining part of this Subsec- 
tion contains slightly modifled details and results of this 
calculation fl. 

Consider the electromagnetic wave of frequency uj = 
c/c, in inflnite array comprised of two types of lossless al- 
ternative a and f3 layers of equal dimensionless thickness 
Aj = 1 with random only dielectric permittivities. Enu- 
merate the layers so that j-th layer occupy the interval 
j — I < z < j and choose all odd layers of a type and 
all even of (3 type. For alternative array, it is natural 
to choose an elementary cell composed of two adjacent 
layers, as the main basic element of the arraji^^^. The 
n-th cell occupies interval 2n — 2 < z < 2n and con- 
sists of (2n — l)-th and 2n layers of type a and /3 cor- 
respondingly. Each layer is characterized by its type a 
magnetic permeability fi^ = 1 = il), refrac- 
tive index ^^^(n) ( i^^(n)), impedance Za{n) = l/vai^) 
{Zf3{n) = ±l/i^Q,(n)), and wave number k^^^ = kua,/3 of 
the wave. 

Within such a model, two systems are considered: the 
H-array when both a and (3 layers are made of right- 
handed materials, and M- array where a layers are right- 
handed material while P layers are of left-handed mate- 
rial. We emphasize that on the contrary of H-stack notion 
where all the layers have the same statistical properties, 
H-array is composed of two different materials with dif- 
ferent statistical properties for odd and even layers. Dis- 
order is incorporated into the model via dielectric per- 
mittivities ea,i3 only, so that refractive index u is a sole 
fluctuation parameter and the upper index in its fluctu- 
ations (^^^^(^) can be omitted 



To calculate Lyapunov exponent of the electromag- 
netic wave of the frequency cj, consider two component 
vector 




with components 



Qn= E{2n-2), Pn = -E'{2n-2) (HI.43) 

UJ 

proportional to the fleld and its derivative at the left edge 
of the n-th cell. These components are real. Therefore 
they automatically correspond to the currentless fleld and 
can be parametrized as 



(HI.44) 



smc 



(compare with Eq. (H.9)). Note that this is currentless 



state in the basis of standing waves while in Section H B 
the basis of running waves was used. 

Using Maxwell equations and appropriate boundary 
conditions on the interfaces of the layers, one obtains 
dynamic equation 



Sn+l = TSn- (III.45) 

Here is the unimodular matrix with elements 



Til = cos (/9q; cos (/:)^ — sin (/9cK sin 

Ti2 = Zq; sin cos (/?^ + cos (/:?CK sin (/:?^, 
T21 = -Z~^ sini^occosip^ - Z^^ cosipocsiiiip^, ' 
T22 = cos (/:)q, cos (/:)/3 — Zc^Z^"*- sin (/pQ, sin 

(HI.46) 

They depend on the cell number n, due to randomized 
refractive indices (HI. 41) entering both the impedances 



Za,i3{'^) and phase shifts cpa^idin)^ 



^«(n) = (^[l + (5^(n)], 
/c^(n) = ±(^[l + (5^(n)]. 



u^in) = 1 + 5a{n), u^in) = ±[1 + 5^{n)]. (HI.41) 

Refractive index fluctuations Sa,f3{n) are assumed to be 
delta-correlated with zero mean value {Sa,i3{n)) = 0, and 
variance a^, 



{6a{n)5^{n')) = a^5a^5nn' , (111.42) 
where angular brackets mean the ensemble average. 



with (p — k/2. 

In On terms, dynamic equations read 



tan 6>n+i 



T2i+r22 tan(9n 
Tii+Ti2tan^^' 



(IH.48) 
(IH.49) 



where now 



23 



1. (Tn + Ti2 tan Of + (Tgi + T22 tan g)^ 
^ ^ 2 l+tan2 6' 



Going to the limit n ^ 00 and using Eqs. (11.24) and 



(III.50) 



(11.25) for localization length = 7 "'we obtain 



(111.51) 



where averaging in the r.h.s. is taken over the stationary 
distribution of the phase 0. 
In the case of weak disorder, 



phase distribution within the main order of perturbation 
theory is uniform 



p{0) = 



(111.56) 



Substit uting this probability density into definition 
(|III.51D and using Eqs. (|III.52D, (|III.53D one gets 



1//^ = 7 = cr^ sin^ (p. 

In the long wave hmit where the phase shift (f is small, 
this result yields the asymptotics 



cr^ < 1 and (acp)'^ <C 1, 



A> 1. 



this distribution p{0) can be explicitly found within the 
framework of a proper perturbation theory. Expanding 
the exact 6>-map ( |III.48 ) up to the second order in pertur- 
batioiP^ and taking into account the uncorrelated nature 
of the disorder (see Eq. (III. 42, one obtains. 



S^{n)U{On 



(III.52) 



This result gives rise to standard A— dependence, ex 
when A ^ oo. In the case of uniform distribution of 5 
over interval [—Qv,Qv\ considered in the Sectio n |IIH it 
exactly coincides with the long wave asymptotic (III. 20) 
of the localization length 

The principally different situation emerges for the M- 
array. In this case the Bloch phase (III. 54) is ide nticall y 
zero. As a result, W{e) = -/7(l9)[r(^)m Eq. (Im^, 



where 



and Eq. (III. 55) leads to a highly nonuniform phase dis- 
tribution 



U{0) = (p -\- sin cp cos(2^ — (/?), 
W{0) = (p[cos{20 - 2(p) ± cos(2i9 - 2(^)] + 
sin (/:)[sin 6 sin(6> - (p) ± sin((9 - (j)/2) sin(6> - (p - 0/2)] + 

sin^ (p sin(46> — 2(p — (j)) cos 0, (III. 53) 

"plus" stands for the H-array, and "minus" refers to the 
M- array, and 



k H-array 
M- array 



(III.54) 



is the unperturbed Bloch phase shift (j) over a unit (a, /3) 
cell 

Now one should write down the Fokker-Plank equation 
related to the dynamic equations (III. 52) 



[U\0)^U\0-m] p{0) 



+2 



dO 



W{0) 



p{0) = 0, 



(III.55) 



find it normalized 7r-periodic solution and calculate aver- 
age in the r.h.s. of (III. 51). 



For H-array, this program can be easily realized. In- 
deed in such a structure the Bloch phase ( HI. 54 ) is non 
zero, and for weak disorder the term in Eq. (HI. 55) con- 



taining (j)/a prevails over the others. Therefore, the 



■sm^ip/U{0). 



(IH.57) 



Fig. p!9| dis plays p e rfect ag reement between analytical ex- 
pressions ( |III.56| ), ( HI. 57) and data obtained by the iter- 
ation of the exact map 



11.45). 



To calculate the Lyapunov exponent via Eq. (HI. 51), 
one needs to perform an ave rage with the distribu- 
tion p{ 0) give n by Eq . ( III.57 K Surprisingly, usage of 

Eqs. ( |IIL51[ ), (|I II.52| ) and ( |III.57| results in zero Lya- 



punov exponent in the main (second order) approxima- 
tion ~ cr^. Therefore, the Lyapunov exponent is deter- 
mined by next orders of the perturbation theory. 

Unfortunately the direct evaluation of high order terms 
in p{0) is rather cumbersome because of huge technical 
complexit}i^. The crucial step which enables authors of 
Ref. [52] to resolve the problem is the following. It is 
known that essential calculation difficulties are often re- 
lated to the non-proper choice of dynamic variables. To 
understand how these variables should be chosen, let us 
analyze the numerical data displayed in Fig. [19) The b- 



panel in this figure demonstrates that the trajectory {i.e. 
the sequence of points (Qn, Pn) has the form of fluctuat- 
ing ellipse specified by angle with respect to axes, and by 
fixed aspect ratio. This results in strongly non-uniform 
phase distribution (d-panel in Fig 19). Therefore, one 
should introduce new variables Qn, Pn by rotating and 
rescaling the axes Q, P, so that the trajectory transforms 
into fluctuating circle. Then, one can expect that the 
distribution of a new phase in the considered approx- 
imation will be uniform. 



24 



RTR- 



(111.58) 



Sn — 



cos Qn 



sin On 



Now the distribution p(0) for new phase c an be f ound 
starting from the quadratic expansion of Eq. (111.48) with 
new coefficients (III. 58) and ^ = 0, 




FIG. 19: (Ref. |52L color online) a) phase space trajectory 
generated by Eq. (|lll.45) for H-array with N — 10"^ , Lp — 



7r/15, for zero disorder (solid circle), and for — 0.003 (scat- 
tered points); b) one trajectory for M- ar ray wit h N = 10^, 
if = 27r/5, (T^ = 0.003 c) p(6 >) from Eq. ( |lll.45| ) for H-array 
(histogram), and Eq. (III. 56) (horizontal line); d) p{0 ) from 
Eq. ( |III.45| ) for M-array (histogram), and Eq. ( |III.57[ ) (sohd 
curve]"! 



©n+l - ©n = [Va{n) - r]^{n)\V{Qn) 
Here the function V{Q) is 



(HL59) 



27? 



V{Q) = sin if) sin(2r — if) sin 20 
I [if — sin if cos(2r — if)] [cos 20 — 1] 

[if + sin if cos(2r - if)] [cos 26 + 1] . (III.60) 



The stationary Fokker-Plank equation corresponding to 
6-map (III. 59) reads 



d 



v\e)^p{e)^v{e)v\e)p{e) 



0. 



From this equation one gets that the phase distribution 
is uniform, p(6) = I/tt, and the trajectory is, indeed, a 
fluctuating circle provided that 



v(e)v'(e) = 0. 



(in.6i) 



With the use of Eqs. ( |III.60 ) and ( |III.61 ), we now can 
obtain the desired expressions for the angle r, parame- 
ter T] and function V{Q) (which is actually no more 0- 
dependent), 



To follow this recipe, let us rotate the vector S ^ S 
RS with the help of unimodular matrix 



R 



COS r sin r 
sin r cos r 



where the angle r describes rotation of the axes in S- 
space, with further rescaling the axes due to free pa ram- 
eter rj. In new coordinates the expressions (HI. 45) and 



(HI. 51) conserve their forms, however, with the rotated 



transfer matrix 



2' 

v{e) 



if -\- sin if 
if — sin if ' 



if"^ — sin^ if. 



(IH.62) 



The results presented in Fig. 20 conflrm success of the 
chosen approach: in new variables the trajectory is a 
fluctuating circle and the phase distribution is uniform. 

T he Lyap unov exponent 7 can be now obtained via 
Eq. (HI. 51) with the change On On- Taking into ac- 



count that 7 vanishes within quadratic appr oximat ion in 
disorder, we expand the 0-map of the form (HI. 48) with 



the coefficients (III. 58) up to the fourth order in per- 



turbation. By substituting the resulting expression into 
Eq. (III. 51 ) and expanding the logarithm within the same 
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Equation ( |IIL63 ) determines the asymptotics for large 
A max(cr, iJJ 



1 

k 



7 



that results in a quite surprising wavelength dependence 
of the localization length, oc X^. Thus, the dependence 
oc A^, numerically found for large A in Refs.^^ and 
confirmed later in should be regarded as the intermediate 
one, apparently emerging due to not sufficiently large 
lengths N over which the average of 7 is performed. 



IV. LOCALIZATION IN COMPLEX MEDIA 

A. Nonreciprocal Transmission in Magnetoactive 
Optical Structures 
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FIG. 20: (Ref. 52 , color online) (a) Phase space trajectory in 
new variables {Q,P); (b) dis tributio n p( Q) gene rated by the 
transformed map with Eqs. ( III. 58 ) and (III. 62), for 7 = 0, 
if 27r/5, cr^ 0.02 and N =10'. 



approximation, after the averaging over with uniform 
distribution, we arrive at final expression 



1 (a^ [{2(f'^ — siii^ (f) cos (p — (p sin (p]"^ 



k 4 

Here the constant 



cp'^ — sin^ cp 



. (III.63) 



2\2 



{Sjnr) - {Sj ny) 

(J2)2 



is specified by the form of distribution of Sa,i3{n). For 
Gaussian and fiat distributions we have ( = 0,-6/5, re- 
spectively. 



In this Subsection we present the results of analyti- 
cal and numerical study of the Anderson localization of 
light propagating through random magnetoactive layered 
structures. We demonstrate that an interplay between 
strong localization and magnetooptical effects produces 
a number of non-reciprocity features in the transmission 
characteristics. 

Magnetooptical effects and nonreciprocity a re wid ely 
exploited in modern optics and applied physic^^^^. In 
particular, magnetoactive periodic structures are cur- 
rently attracting growing attention^^ The main phe- 
nomena of interest are the enhanced Faraday effect on 
resonances^ and one-way propagation (nonreciprocal 

'86ti89l 



transmissionp^^*^ employed for the concept of optical 
insulators. The resonant Faraday effect has also been 
shown in connection with the localization of light in ran- 
dom layered structure^^. 

Here we examine the transmission properties of one- 
dimensional random layered structures with magneto- 
optical materials. We employ short-wavelength approx- 
imation, where the localization is strong, and consider 
both Faraday and Voigt geometries. In the Faraday 
geometry, magneto-optical correction to the localiza- 
tion length / results to a significant broadband non- 
reciprocity and polarization selectivity in the typical, ex- 
ponentially small transmission. In the Voigt geometry, 
averaging over random phases suppresses the magneto- 
optical effect, in contrast to the case of periodic struc- 
tures where it can be quite pronounced^^ At the 
same time, in both the geometries we reveal the non- 
reciprocal frequency shifts of narrow transmission reso- 
nances, corresponding to the excited localized states in- 
side the structure^ This offers efficient unidirec- 
tional propagation at the given resonant frequency. 

Consider the light transmission through the long stack 
composed with magnetooptical materials in the short- 
wavelength approximations. In the localized regime, we 
can neglect in Eq. (11.21) the external interface transfer 



matrices F^", F^^ just replacing the exact matrix f by 
the truncated matrix T' 
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T' = FnSnFn-iSn-iFn-2 ••• F2S2F1S1. (IV.l) 

Then, if the wavelength within the k-th layer is much 
shorter than the variance of the layer thicknes^^, then 
the ph ases (fk niodulo 27r in the propagation matrices Sj 
(11.20) are independent and nearly uniformly distributed 



in the range (0, 27r). In this approximation, the transmit- 
tance corresponding to the transfer matrix (IV.l) after 



averaging over all phases (fk, is reduced to the product of 
the transmittances of separate layers^^ and, furthermore, 
to the product of transmittances of the interfaces only^ 



2N 



ln(r)^^lnr, , r, = (IV.2) 



Substitution of Eq. ( |IV.2| ) into Eq. ( |II.1[ ) in the limit 
N ^ 00 yields the simple expression for the localization 
length 



in the short-wavelength approximation. 

This result can be easily extended to any number of 
alternating layers. For instance, considering a random 
structure consisting of three types of alternating layers, 
'a', and '7', one has 



^Pa(3^ ^pf^l^ ^p70i^ 



Transmission through a one-dimensional lossless linear 
medium is always reciprocal if there is only one (but prop- 
agating in two directions) mode in the system. Indeed, 
while the forward transmission of the wave incident from 
the left on the medium is described by the 2x2 trans- 
fer matrix T with transmission coefficient T and trans- 
mittance T, the backward transmission of the reciprocal 
wave incident from the right is characterized by the in- 
verse transfer matrix with the same transmission 
coefficient and transmittanc^^^. 

If the system possesses two or more uncoupled modes 
labeled by index the waves are marked by the propaga- 
tion direction v and mode indices: h'^'^ . Still, the forward 
and backward propagation of each mode ^ through the 
system with incident waves of types (+, ^) and (— , ^) are 
described by the 2x2 transfer matrices and (T^)~^ 
characterized by the same transmittance . However, 
the wave reciprocal to (+, is determined by the time- 
reversal operation which changes v —v (because of the 
k ^ — k transformation) but can also affect In par- 
ticular, if the time reversal operation changes the sign of 
the mode index: ^ ^ —^^ then the reciprocal wave will be 
( — , — (^) rather than the backward wave of the same mode, 
(— ,c^). Accordingly, the transmittances of the mutually 
reciprocal waves through the system, and ^ can 
be different. This signals nonreciprocity in the system. 




etizatio^ 




FIG. 21: (Ref. 34 , color online.) Schematic picture of the 
wave transmission and reflection from a random-layered struc- 
ture consisting of two types of alternating layers 'a' (here 
- magnetoactive material) and (here - air) with random 
widths. Magnetization of the medium, wave polarizations and 
directions of propagation are shown for the Faraday and Voigt 
geometries. 



Non-reciprocity in the system under consideration 
originates from the difference between the modes <; and 
— and does not depend explicitly on the direction of 
incidence v. Therefore, in practice, it is sufficient to 
compare only forward transmissions of the modes de- 
scribed by the transfer matrices T^"^ and transmittances 

There are two main geometries typical for magneto- 
optical problems^^: the Faraday geometry, where the 
magnetization is collinear with the direction of propa- 
gation of the wave, and the Voigt (or Cotton-Mouton) 
geometry, where the magnetization is orthogonal to the 
direction of propagation of the wave (see Fig|2T]). Below 
we study the averaged transmission decrement and in- 
dividual transmission resonances in both geometries and 
show that propagation of light in disordered magnetoac- 
tive layered media offers nonreciprocal transmission. 

In the Faraday geometry both magnetization and the 
wave vector are directed across the layers, i.e., along the 
z-axis (see Fig 21). We assume that the magnetic tensor 



is equal to one and the magneto-optical effects are de- 
scribed exclusively by the dielectric tensor which in the 
Faraday geometry has the form^^ 



£ 
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The eigenmodes of the problem are circularly polarized 
waves of magnetic H 
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V2 

and electric E 




^i{vkz—ujt) 



^;,c^ = ±l, (IV.4) 



k 



(IV.5) 



fields. Here H^'"^ {E'^'^) are the wave amplitudes, whereas 
k is the propagation constant affected by the magnetiza- 
tion parameter q and depending on 



LO Q 

C £ 



ko 



(IV.6) 



In the linear approximation in g', /c :^ nko{l -\- ^q/2). 

Parameter ^ is the mode index which determines the 
direction of rotation of the wave field. In this manner, 
the product represent the helicity 



X 

which distinguishes the right-handed (x = +1) and left- 
handed (x = —1) circular polarizations defined with 
respect to the direction of propagation of the wave. 
Note that the time reversal operation keeps helicity un- 
changed, whereas ^ changes its sign^^. Thus, the recipro- 
cal wave is given by H~^'~'' , precisely as described above. 

The total fi eld in a layer is the sum +H~''' of the 
eigenvectors (IV.4) with the amplitudes H^'"^ . Consider 



the wave transformation at the interface between the me- 
dia 'a' and 'h\ The helicity of the wave flips upon the 
reflection and remains unchanged upon transmission. As 
a result, parameter ^ remains unchanged, so that there is 
no coupling between the modes with <^ = +1 and ^ = —1 
(see Fig 21), and these modes can be studied indepen- 



dently. From now on, for the sake of simplicity, we omit 
in superscripts and write explicitly only the values of 

the direction parameter v = ^1. 

Using the standard boundary conditions for the wave 

electric and magnetic fields at the interface for the 

normalized fields 



H- 



L, = F^^ hs 



with the normalized interface transfer matrix 



2y/kak(3 



k(3 — key kf^ + ko^ 



(IV.7) 



where kc^^^ are the wave numbers (IV.6) in the corre- 
sponding media. 
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FIG. 22: (Ref. [34 , color online) Localization decrement Hi 
vs. magneto-optical parameter Q for opposite modes prop- 
agating through a two-component random structure in the 
Faraday geometry (see details in the text). The modes with 
^ = zbl correspond to either opposite circular polarizations or 
propagation directions. Numerical simulations of exact equa- 
tions (symbols) and theoretical formula (IV.8) (lines). 



Calculating the localization decrement from Eq. (IV. 3 ) 



with Eqs. (IV.6) and (IV.7), we obtain in the linear ap- 



proximation in q\ 
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Thus, the localization decrement acquires the first-order 
magneto-optical correction n^^^ caused by the Faraday 
effect. This correction depends on i.e., on the polariza- 
tion helicity x and the propagation direction v through 
= X^- For the reciprocal waves with the same x and 
opposite n^^^ has opposite signs. This signals nonre- 
ciprocal localization in a Faraday random medium. In 
practice, the nonreciprocal difference in the transmission 
decrements (IV.8) can be observed by changing sign of 



either propagation direction v (with the helicity fixed), 
or polarization x, or magnetization q. 

Despite the magneto-optical correction to the localiza- 
tion decrement is small in magnitude, k!^^^ <C ti^^\ it 
still might result in a significant difference in the typi- 
cal transmission spectrum. This difference is described 
by an additional factor of ex exp[— 2A^a>:^^^] in transmit- 
tance, which is exponential with respect t o the length 
of the structure. Hence, small correction (IV.8) brings 



about significant broadband nonreciprocity or polariza- 
tion selectivity in the typical small transmission when 

Fig. 22 shows dependence of the localization decre- 
ment on the magnetization parameter Q = eq calculated 
numerically and compared to analytical result (IV.8). 



Numerical simulations were performed for the struc- 
ture containing M = 2N = 90 alternating layers of 
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air {e = 1,(5 = 0), and bismuth iron garnet (BIG), 
with dielectric constant e = 6.25 and magneto-optic 
parameter reaching Q = 0.06. The thicknesses of layers 
were randomly distributed in the range 50^150 /im 
(i.e., w = 100/im, d = 50/im), whereas the excitation 
wavelength was 632 nm. The averaging was performed 
over 10^ realizations of the random sample. One can see 
excellent agreement between numerical simulations and 
analytical results showing linear splitting of the = 1 
and = — 1 localization decrements as a function of the 
magneto-optic parameter. 

In the Voigt geometry, the dielectric tensor i^ 
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The first-order interaction of the wave with the magneti- 
zation occurs only upon oblique propagation of the wave 
in the xz-plane, i.e., when kx = const ^ (see Fig 21). 

The eigenmodes of the problem are the TE mode which 
is uncoupled from the magnetization, and TM mode with 
the tangential components 



Here parameters w = ±1 and <r = ±1 indicate prop- 
agation in the posit ive and negative z and x directions, 
respectively, k\\ = i/fc^ — fc^, k±_ = \kx\, whereas 



ic;qk±_ + vk\\ 
e{l-q^)ko 



(IV.9) 



k — nko\/l — cp' 



In the linear approximation in q^ A^''^ (^;/c|| + 
i^qk±)/{£ko) and k nko so that the magnetization 
affects imaginary parts (i.e., phases) of the amplitudes 
A^''' and d oes n ot aff ect the propagation constant, cf. 
Eqs. drvTsl) and (IrvTel). 



In the Voigt geometry, direction of the transverse wave 
vector component, serves as the mode index. The mu- 
tually reciprocal waves are H^'^ and H~^'~^ because the 
time reversal transformation reverts the whole wave vec- 
tor, k — k. 

The parameter ^ is not changed upon reflection and 
transmission through the layers, i.e., the modes with 

= ±1 are uncoupled from each other. Therefore, for 
the sake of simplicity, we omit the mode index in super- 
scripts, and write explicitly only the values of the direc- 
tion parameter = ±1. 

Using the standard boundary conditions for the wave 
electric and magnetic fields at the interface, for the 

normalized interface transfer matrix we obtair 



FIG. 23: (Ref. [34 , color online.) Transmission spectra of 
a random magneto-optical sample in the Faraday geometry 
(see details in the text) for waves with = ±1. While the 
averaged localization decrements are only slightly different 
(Fig |22[ ), all individual resonances are shifted significantly as 
compared with their widths, Eq. (IV. 10). 
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In contrast to the Faraday geometry, in the Voigt 
geometry the linear magneto-optical correction changes 
only phases of the transmission and reflection coeffi- 
cients, whereas corrections to the interface transmit- 
tance start with the terms (X q^. In short-wave limit, 
only these transmittances determine the total transmit- 
tance, Eq. (IV. 2). Therefore, a short-wavelength trans- 



mission through a random multilayered stack is recip- 
rocal and is not affected by magnetization in the first- 
order approximation. In the short-wave limit, this state- 
ment remains true for any number of types of alternat- 
ing layers. It was verified numerically for the three-layer 
systenP^. At the same time, a periodic structure with a 
cell consisting of three different layers (which breaks the 
mirror reflecti on sy mmetry) can demonstrate significant 
nonreciprocit}ff^I^ but beyond the short-wave approxi- 
mation. 

Averaged localization decrement is associated with ex- 
ponential deca y of t he incident wave deep into the infi- 
nite sample^ 49l56] ^ finite sample, this is so only 
for typical realizations. However, there exist some res- 
onant realizations of the sample at a given frequency 
(or, equivalently, resonant frequencies for a given real- 
ization) where transmission is anomalously high and is 
accompanied by the accumulation of energy inside the 
sample. Such resonant transmission corresponds to 
excitation of the Anderson localized states (quasi-modes) 
inside the sample. Akin to the resonant localized states 
in photonic crystal cavities, the transmission resonances 
in random structures are extremely sensitive to small per- 
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turbations: realization^, absorptionp^ nonlinearity,'^ 
and, as we show here, magnetoactivity. 



of its maximum and minimum correspond to the width 
of the original resonance and are practically unchanged. 
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FIG. 24: (Ref. [3T, color online.) Differential transmittance, 
— 7 , for two resonances from Fig 22 as dependent on the 
value of magneto-optical parameter Q, cf. Eq. (IV. 10). 



Figure [23] shows transmission spectra for two modes 
= ±1 (i.e., either with opposite helicities or propaga- 
tion directions) in one realization of a magnetooptical 
sample in the Faraday geometry. The parameters of the 
sample are the same as in Section IIIA with Q = 0.06. 
One can see strong splitting of the = ±1 transmission 
resonances which have exponentially narrow widthJ^ 
(X A>:exp(— /^7V)/2'iD . This offers strongly nonreciprocal, 
practically unidirectional, propagation or polarization se- 
lectivity in the vicinity of resonant frequencies. 

To estimate the splitting of resonances, we note that 
the wavenumbers in magnetooptical materials are shifted 
due to the Faraday effect, Eq. (IV.6). Hence, the shifts 



of the resonant wavenumbers of the random Faraday 
medium can be estimated by averaging of this shift over 
different materials in the structure: 



qnko 



(IV. 10) 



where (...) stands for some average of (...). Using qn ^ 
{Qa^a + Qb^b)/"^ for estimation in the two-component 
structure, we obtain AAres ^ — ^^3.6 nm, which agrees 



with the ^-dependent splitting observed in Fig. 23 



Figure [24| displays the differential transmission for the 
waves with = +1 and lying in a narrow frequency 
range in Fig. 23 In agreement with estimation (IV. 10), 



one observes the linear dependence of the resonance split- 
ting on magnetization. 

In the Voigt geometry, the resonances also allow non- 
reciprocal transmission and demonstrate splitting of the 
resonant frequencies. In Fig. 25, the differential trans- 
mission is shown for reciprocal waves with <^ = ±1 in 
the vicinity of one resonance for the three-component 
structure considered in Section IIIB. The splitting is 
very small in this case, and <^ = +1 and = —1 res- 
onances overlap significantly. Because of this, the dif- 
ferential transmittance in Fig. 25 is tiny, its amplitude 
linearly grows with Q, whereas the frequency positions 
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FIG. 25: (Ref. [34 , color online.) Differential transmittance, 
— 7 , for in the vicinity of a single resonance in the Voigt 
geometry (see Section IIIB for details) as dependent on the 
magneto-optical parameter Q. 

Unlike the wave- number shift in the Faraday geom- 
etry, the nonreciprocal shift of resonant frequencies in 
the Voight geometry arises from the phases of the ampli- 
tudes A, Eq. (IV.9). These phases are responsible for the 



phases of transmission coefficients between the layers and 
can be estimated as (/) ~ q{(;k^) / {vk\\) = gtan^, where 
is the angle of propagation with respect to the z-axis. 
The phases accumulated at a layer effectively shift the 
wave numbers as vAk\\ = A/ccos6> ~ (l>/w^ where w is 
the thickness of the layer. Averaging over different ma- 
terials in the random layered structure, we estimate the 
nonreciprocal shift of the resonant wave number: 



Akr, 



q smt 



q\ smt 



w cos^ w cos^ 
This shift is (^-dependent, i.e., nonreciprocal, and much 
smaller than the Faraday-geometry shift (IV. 10) as kw > 
kd ^ 27T in the short- wavelength limit. For the pa- 
rameters in use, with Q = 0.06, we have AAres ^ 



-(^3 • 10- 
Fig. 25. 



nm, which agrees with the data plotted in 



B. Charge Transport in Disordered Graphene 

Shortly after the discovery of highly unusual physi- 
cal properties of graphene it was realized that the elec- 
tron transport in this material had many common fea- 
tures with the propagation of light in dielectrics. His- 
torically, the analogy between Maxwell equations and 
those used in the relativistic electron theory has been 
discussed in different contexts and for various purposes 
(see, for example, S^'^^M^ since 1907 when Maxwell equa- 
tions were reduced^ to an alternative, more concise form 
by introducing a complex field F = E + zH: 



(IV.ll) 
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where ^ is the 4- vector with components —Fx -\-iFy^ 
Fz, Fx ^ iFy^ n is the refraction index, and the com- 
ponents of the 3- vector p, are the Dirac matrices p^, 
i = 1,2,3 (Pauh matrices in which the units are replaced 
by the unit 2x2 matrices). 

In the last few yeas, this activity perked up due to the 
recent developments in the physics of graphene. Nowa- 
days it is well understood that under some (rather gen- 
eral) conditions, Dirac equations describing the charge 
transport in a graphene superlattice created by applying 
inhomogeneous external electric potential could be re- 
duced to Maxwell equations for the propagation of light 
in a dielectric medium. To better understand the physics 
of charge transport in graphene subject to a coordinate- 
dependent potential, in what follows, we compare the re- 
sults for graphene with those for the propagation of light 
in layered dielectric media (for more analogies between 
quantum and optical systems, see, e.g. Ref.^^^^). Addi- 
tional analogies, not discussed here, also exist with the 
transport and localization of phonons in different kinds 
of periodic and random one-dimensional structure^^^Hinni^ 

As it was shown above, the light transport of elec- 
tromagnetic waves in multilayered media is described in 
terms of the transfer matrices of two types. The first type 
is formed by diagonal matrices Sj corresponding to the 
propagation of wave throug h the j -th layer. These matri- 
ces are the same as in Eq. (11.20) (up to the signs of the 
exponent). The second type is formed by the interface 
transfer matrices describing transformation of the 

amplitudes of the electromagnetic waves at the interface 
between j— th and {j + 1)— th layers and having the form 



1 



2 cos Oj-^i 



(IV.12) 



where 



Gjj^i = cos^,+i ± cos^, • sgn(^,^,+i)^ (IV. 13) 
for s-polarized waves and 



cos^j+i ± cos^j • sgn(^'j^'j+i) (IV. 14) 



for p-polarized waves. Here, Oj is the angle of the 
propagation within the j— th layer, Zj and Vj are the 
impedance a nd the refractive index of the jth layer de- 
fined by Eq. (III.l). Signs ± correspond, respectively, to 
R— and L— dielectric layers with positive and negative 
refractive indices. 

In the case of the charge transport in a graphene su- 
perlattice created by a piecewise-constant electrostatic 
potential depending on one coordinate x in the plane 
(x, y) of the graphene layer, the analogues transfer ma- 
trix, which describes the transition through the interface 
between adjacent regions with different values of the po- 
tential, has the form^^ 



1 



2 cos e 



J+1 



t(-) 



(IV.15) 



where 



■^gn[{e-Uj){£-Ujj^i)], (IV.16) 



Here the 



is given by equation tan^j 



where (3 is the projection of the 



dimensionless momentum on y axis, e and Uj are the 
dimensionless energy of the charge carrier and the scalar 
potential of the j-th layer. If Oj is real, it coincides with 
the angle of the propagation of electron within the j-th 

layer. 

Comparison of Eqs. (IV.12) and (IV.15) shows that 



the role of the refractive index v in graphene is played 
by the difference e — u. In particular, a layer, in which 
the potential exceeds the energy of the particle, u> e/\s 
similar to a L-slab with negative refractive index (meta- 
material), while a layer where u < e/\s similar to normal 
material. It is due to this similarity that a junction of 
two regions having opposite signs oi u — £ (so-called p-n 
junction) focuses Dirac electrons in graphen^^, in the 
same way as an interface between left- and right-handed 
dielectrics focuses electromagnetic waves^^ . 

This analogy is not complete: although the equations 
are akin, the boundary condition s are, generally speak- 
ing, different. As a result Eq. ( |IV.16| ) (for graphene) 
does not contain factor Zjj^i/Zj which is present in 



Eqs. (IV. 13), (IV. 14) and determines the reflection co- 
efficients at the boundary between two dielectricJ^^. 
Another important distinction between transfer matri- 
ces JF (graphene) and F (electromagnetic waves) is that 

is a complex-valued matrix, while the F is always 
real. This is manifestation of the fundamental differ- 
ence between graphene wave functions and electromag- 
netic fields in dielectrics. The graphene wave functions 
are complex-valued spinors which describe two different 
physical objects: particles (electrons) and antiparticles 
(holes). The electromagnetic fields are real that reflects 
the fact that photons do not have antiparticles (antipho- 
ton is identical with photon). These distinctions bring 
about rather peculiar dissimilarities between the conduc- 
tivity of graphene and the transparency of dielectrics. 

However in the particular case of normal incidence 
9j = Ojj^i = and equal impedances Zj = ^j+i, the 
transmission of Dirac electrons through a junction is 
similar to the transmission of light via an interface be- 
tween two media with different refractive indices (but 
equal impedances) . Such an interface is absolutely trans- 
parent to light and therefore both p-n and p-p junc- 
tions are absolutely transparent to the Dirac electrons in 
graphene^^"-*^^^". This is related to the absence of backscat- 
tering and antilocalization of massless Dirac fermions 
caused by their spin propertie d^^ l ^^l This also explains 
Klein paradox^^^ (perfect transmission through a high 
potential barrier) in graphene systems, and leads (to- 
gether with symmetry and spectral flow arguments) to 
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the surprising conclusion that Dirac electrons are de- 
localized in disordered ID graphene structure, provid- 
ing a minimal non-zero overall transmission, which can- 
not be destroyed by fluctuations, no matter how strong 
they are^^^. However, this statement (being correct in 
some sense) should be perceived with a certain caution. 
Indeed, many features of Anderson localization can be 
found in random graphene systems. It has been shown 
that although the wave functions of normally inci- 
dent {0 = 0) particles are extended and belong to the 
continuous part of the spectrum, away from some vicin- 
ity of ^ = 0, 1-D random graphene systems manifest all 
features of disorder-induced strong localization. In par- 
ticular, for a long enough, disordered graphene superlat- 
tice the transmission coefficient, T, as a function of the 
angle of incidence, ^, (or of the energy ifO ^ Ois fixed) 
has typical for Anderson localization shape. Fig. l26| 




T(0) 




FIG. 26: (Ref. 47 ) Transmission coefficient T(0) for periodic 
(thin black line) and disordered (bold blue line) graphene. 

Along with continuous of typical angles (or energies), 
for which the transmission is exponentially small, there 
exists a quasi-discrete random set of directions where the 
sample is well transparent, i.e., the transmission coeffi- 
cient is close to one. At these angles, the wave functions 
are exponentially localized (Fig. 27), with the Lyapunov 
exponent (inverse localization length l^) being propor- 
tional to the strength of disorder. 

Charge transport in a graphene sheet subjected to a 
disordered electrostatic potential is determined by the 
ratio between its values u{^) and the energy e of the 
particle. In particular, in randomly-layered potential 
Uj = uo{j) + Auj {j is the number of a layer, uo{j) 
is a non-random function, Auj are independent ran- 
dom variables homogeneously distributed in the interval 
[—Su, 6u]) it manifests essentially different features in the 
following three different systems^^: 



(i) 



< £, 1^0 (j) is a periodic function. In this 



case, a relatively weak disorder drastically changes the 



FIG. 27: (Ref. ^47|) Spatial distribution of the wave func- 
tion localized inside the sample for 6 marked by red arrow in 
Fig. [26] 



transmission spectrum. All features of the spectrum of 
the underlying periodic structure are washed out, and a 
rather dense (quasi-) discrete angular spectrum appears, 
with the corresponding wave functions being localized at 
random points inside the sample (disorder-induced res- 
onances). However, there is one fundamental difference 
from the usual Anderson localization: in the vicinity of 
normal incidence, the transmission spectrum of graphene 
is continuous with extended wave functions, and the 
transmission coefficient is finite (T = 1 at 6> = 0). It is 
this range of angles that provides the finite minimal con- 
ductivity, which is proportional to the integral of T{6) 
over all angles 0. 

(ii) e < uo{j) = const. Under these conditions, the 
transmission of the unperturbed system is exponentially 
small (tunneling) and, rather unusually, gets enhanced 
by the fluctuation of the potential. 

(iii) £ = 0, uo{j) is a periodic set of numbers with 
alternating signs. The behavior of the charge carriers in 
the potential of this type is most unusual. It is character- 
istic of two-dimensional Fermions and have no analogies 
in electron and light transport. The disorder obliterates 
the transmission peaks of the underlying periodic system, 
makes much wider the transparency zone around normal 
angle of incidence, and gives rise to a new narrow peak 
in the transmission coefficient, associated with wave lo- 
calization in the random potential. Unlike the peaks in 
the periodic structure, the wave function of this disorder- 
induced resonance is exponentially localized. In distinc- 
tion to the case (i), the transmission in (iii) is extremely 
sensitive to fluctuations of the applied potential: relative 
fluctuations Au/uq = 0.05 reduce the angular width of 
the transmission spectrum more than four times. 

Propagation of light in analogous L-R and R-R disor- 
dered dielectric structures demonstrates completely dif- 
ferent behavior. As the degree of disorder (variations of 
the refractive index) grows, the averaged angular spec- 
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tra quickly reach their asymptotic "rectangular" shape: 
a constant transmission in the region where all interfaces 
between layers are transparent followed by an abrupt de- 
crease in transmission in the region of angles where the 
total internal reflection appears. 



C. Bistability of Anderson Localized States in 
Nonlinear Media 

Recent renewed interest to Anderson localization is 
drive n by a series of experimental demonstrations in op- 
tics^^H^and Bose-Einstein condensates^^ One of the 
important issues risen in these studies is that the disor- 
dered systems can be inherently nonlinear, so that an 
intriguing interplay of nonlinearity and disorder could be 
studied experimentally. 

Nonlinear interaction between the propagating waves 
and disorder can significantly change the inter ference ef- 
fects, thus fundamentally affecting localizatio ri^^ * ^^ * ^ -^^^. 
However, most of the studies of the localization in ran- 
dom nonlinear media deal with the ensemble- averaged 
characteristics of the field, such as the mean field and 
intensity, correlation functions, etc. These quantities de- 
scribe the averaged, typical behavior of the field, but they 
do not contain information about individual localized 
modes (resonances), which exist in th e localized re gime in 
each realization of the random samplj^ ^^ l ^^ l ^-*^ ! -*^^ -^. These 
modes are randomly located in both real space and fre- 
quency domain and are associated with the exponential 
concentration of energy and resonant tunneling. In con- 
trast to regular resonant cavities, the Anderson modes 
occur in a statistically-homogeneous media because of the 
interference of the multiply scattered random fields. Al- 
though the disorder-induced resonances in linear random 
samples have been the subject of studies for decades, the 
resonance properties of nonlinear disordered media have 
not been explored so far. 

In this Section we present the study of the effect of 
nonlinearity on the Anderson localized states in a one- 
dimensional random medium^^. As a result of interplay 
of nonlinearity and disorder, the bistability and nonre- 
ciprocity appear upon resonant wave tunneling and ex- 
citation of disorder-induced localized modes in a manner 
similar to that for regular cavity modes. At the same 
time, weak nonlinearity has practically no effect on the 
averaged localization background. 

First, let us consider a stationary problem of the 
transmission of a monochromatic wave through a one- 
dimensional random medium with Kerr nonlinearity. 
The problem is described by the equation 



(IV. 17) 



where ip is wave field, x is coordinate, k is wave number 
in the vacuum, n = n{x) is the refractive index of the 
medium, and x is the Kerr coefficient. 

In the linear regime, xlV^P = 0? the multiple scat- 
tering of the wave on the random inhomogeneity n'^{x) 



brings about Anderson localization. The main signa- 
ture of the localization is an exponential decay of the 
wave intensity, / = deep into the sample and, 

thus, an exponentially small transmissior l^ l ^ l ^ *^ ^olt'^ ^ 
/in exp(— 21///) <^ 1. Here L is the length of the sample 
and / is the localization length which is the only spatial 
scale of Anderson localization. Along with the typical 
wave transmission, there is an anomalous, resonant trans- 
mission, which accompanies excitation of the Anderson 
localized states inside the sample and occurs at random 
resonant wave numbers k = k^eso^ - In this 

case, the intensity distribution in the sample is charac- 
terized by an exponentially localized high-intensity peak 
inside the sample, /peak ^ /in, and a transmittance much 

higher than the typical one: I^^^^ ^ ^olrt^^- 

Excitation of each localized mode inside the random 
sample can be associated with an effective resonator cav- 
ity located in the area of field localization and bounded 
by two potential barriers with exponentially small trans- 
parencies^. According to this model, the transmittance 
spectrum T{k^Iout) in the vicinity of a resonant wave- 
length for the c ase of w eak nonlinearity (xlV^P <^ 1) is 
given in the forn j-^^ * ^-^ ^ 



T(fe,/out) 



hn 1 + [Ax/out + ^] 



2 ' 



(IV.18) 



where Tres is the transmission coefficient at resonance, 
and dimensionless parameters A and 5 characterize, re- 
spectively, the strength of the nonlinear feedback and the 
detuning from the resonant wave number: 
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. (IV.19) 



Equation ( |IV.18 ) establishes relation between the in- 
put and output wave intensities, which is given by a 
cubic equation with respect to /out- It has a universal 
form typical for nonlinear resonators possessing optical 
bistabilit^^^^. From Eq. ( |IV 19| ) it follows that in the 
region of parameters: 



AJ<0, J2>3, |x|/in> 



1 



3^3 \A\T,, 



the dependence /out(^in) is of the S-type and the sta- 
tionary transmission spectrum T(k) is a three-valued 
function. In most cases, one of the solutions is unsta- 
ble, whereas the other two form a hysteresis loop in the 
^out(^in) dependence (see Figs. 28 and 29). 

It is important t o emph asize two features of the equa- 
tions ( |IV.18 ) and ( |IV.19 ), describing the nonlinear res- 
onant transmission through a localized state. First, 
they have been derived without any approximations 
apart from the natural smallness of the nonlinearity and 
Lorentzian shape of the spectral line. Second, although 
the resonant transmission, the effect of nonlinearity, and 
bistability owe their origin to the excitation of the Ander- 
son localized mode inside the sample, equations (IV.18) 
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FIG. 28: (Ref. 31 , color online.) Nonlinear deformations of 
the transmission spectra of two random resonances at differ- 
ent intensities of the incident wave. Numerical simulations 
of the Eq. 



IV.17 



IV.18 



(sym- 



(curves) and theoretical Eq. 
bols) are shown for the case of defocusing nonlinearity, x > 0- 
Light-grey stripes indicate three-valued regions for the high- 
intensity curves, where only two of them (corresponding to 
the lower and upper branches) are stable. 



and (IV. 19) contain only quantities which can be found 
via outside measurements^^. 



Figure 28 shows nonlinear deformations of the reso- 
nant transmission spectra T{k) for different values of /in, 
which exhibit transitions to bistability. The analytical 
dependence T{k) given by Eqs. (IV.18, IV. 19) with the 
parameters T^es, Q, and A found from the numerical ex- 
periments are in excellent agreement with the direct nu- 
merical solutions of Eq. (IV. 17)^^^. In numerical simu- 



lations of stationary regime we used the standard 4-th 
order Runge-Kutta method. We note, that the incident 
field amplitude is a single-valued function of the trans- 
mitted field. Thus, we solve second-order ordinary differ- 
ential equation Eq. (IV.17) using transmitted field value 



as the boundary conditions for the equation. 

The dimensionless parameters Tres, and Q from Eqs. 
(IV.18, IV. 19), can also be estimated from a simple res- 



onator model of the Anderson localized states 



9 78 91. 



4T1T2 
(T1+T2)' 



T1+T2 
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where 



Ti exp [-2Xres/l] , T2 exp [-2(L - Xres)/l] 

are the transmission coefficients of the two barriers that 
form the effective resonator, x^es is the coordinate of the 
center of the area of field localization, / is the localization 
length, and L is the length of the sample. 

Introducing a weak Kerr nonlinearity into the res- 
onator model, one can also estimate the nonlinear feed- 
back parameter A: 



A - Q/n^T2 



(IV.21) 



where is the mean value of n^(x). 

It is important to note that each disorder-induced res- 
onance is associated with its own effective cavity, so that 




2 4 6 

Incident power, xli^{x10'^) 



2 3 4 

Time (10 -""^s) 



FIG. 29: (Ref. color online.) Stationary and FDTD 

simulations showing hysteresis loops in the output vs. input 
power dependence for three different resonances. Panel (d) 
shows deformation of the transmitted Gaussian pulse corre- 
sponding to the hysteresis switching on the resonance 2. 



the disordered sample can be considered as a chain of 
randomly loc ated coupled r esonators^ ^. 

Equations (IV.20, IV.21| ) enable one to estimate the 



values of the parameters describing the no nlinear reso- 
nant wave tunneling in Eqs. 
only the basic parameters of t : 
ization coordinate and th e localization length. Inpar- 
ticular, substituting Eqs. (|IV.20l |IV.21|) into Eq. (11 V. 191) 



IV.20, 


IV.21 ) by knowing 


le loca 


ization - the local- 



and taking into account that the most pronounced trans- 
mission peaks correspond to the localized states with 
X L/2 and Ti ~ T2, we estimate the incident power 
needed for bistability of localized states 



lx|/ir 



> 



exp(-2L//) 



■^resO 



(IV.22) 



For the parameters used in our simulations this gives 
quite reasonable value |x|^in ^ 10~^. If we increase 
the length of the sample, the Q-factors of the resonances 
grow, and the incident power needed to observe the bista- 
bility becomes smaller. 

To demonstrate temporal dynamics upon the bistable 
resonant tunneling, an explicit iterative nonlinear finite- 
difference time-domain (FDTD) scheme was imple- 
mented. For precise modeling of the spectra of narrow 
high-Q resonances, fourth-order accurate algorithm was 
used, both in space and in time^^^, as well as the Mur 
boundary conditions to simulate open boundaries and 
total- field/scattered- field technique for exciting the inci- 
dent wave. Sufficient accuracy was achieved by creat- 
ing a dense spatial mesh of 300 points per wavelength 
{dx = A/300). To assure stability of the method in non- 
linear regime, the time step was selected as dt = dx/3c^ 
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FIG. 30: (Ref. [31], color online) (a) Non-reciprocal trans- 
mission through the nonlinear disordered structure, showing 
different output powers for identical waves incident from dif- 
ferent directions, (b) Corresponding shape of the incident 
pulse, and pulses transmitted in different directions. 



and each simulation ran for = 2 * 10^ time steps. To 
compare the results of the FDTD simulations with the 
steady-state theory, the transmission of long Gaussian 
pulses with cen tral fre quencies and amplitudes satisfy- 
ing conditions (IV. 22) was considered, see Fig. 29 ^d). 
With an appropriate choice of the signal frequencies, 
we observe hysteresis loops in the /out(^in) dependencies 
which are in excellent agreement with stationary calcu- 

Characteristic transi- 



lations, as shown in Figs. 29 



a-c 



tional oscillations accompany jumps between two stable 
branches, and strong reshaping of the transmitted pulse 
evidences switching between the two regimes of trans- 
mission. Fig. 29 'd). We note, that different choice of the 
signal frequencies near the resonance can lead to various 
other behaviors of output vs. input curves, with trans- 
mission either increasing, when nonlinear resonance fre- 
quency shifts towards the signal frequency, or decreasing 
in the opposite case. 

In addition to the bistability, the resonant wave tun- 
neling through a nonlinear disordered structure is nonre- 
ciprocical. As is known for regular systems, nonsymmet- 
ric nonlinear systems may possess nonreciprocal trans- 
mission properties, resembling the operation of a diode. 
An all-optical diode is a device that allows unidirectional 
propagation of a signal at a given wavelength, which may 
become useful for many applicationJ^^. A disordered 
structure is naturally asymmetric in the generic case, and 
one may expect a nonreciprocal resonant transmission in 
the nonlinear case. To demonstrate this, we modeled 
propagation of an electromagnetic pulse impinging the 
same sample from different sides and monitor the trans- 
mission characteristics. One case of such nonreciproci- 
cal resonant transmission is shown in the panel (a) of 
Figj30| We observe considerably different transmission 



properties in opposite directions with the maximal in- 
tensity contrast between two directions 7.5:1. Moreover, 
the threshold of the bistability is also significantly differ- 
ent for two directions: there is a range of incident powers, 
for which the wave incident from one side of the sample 
is bistable, while there is no signs of bistability for the 
incidence from the other side. Figure 30 ^b) shows the 
pulse reshaping for incidence from opposite sides of the 
structure. 



In this Section we have presented the study of new 
manifestations of the interplay between nonlinearity and 
disorder. It is shown that even weak nonlinearity affect 
dramatically the resonant transmission associated with 
the excitation of the Anderson localized states leading to 
bistability and nonreciprocity. Despite random character 
of the appearance of Anderson modes, their behavior and 
evolution are rather deterministic, and, therefore, these 
modes can be used for efficient control of light similar to 
regular cavity modes. These results demonstrate that, 
unlike infinite systems, the Anderson localization in finite 
samples is not destroyed by weak nonlinearity - instead 
it exhibits new intriguing features typical for resonant 
nonlinear systems. 



V. CONCLUSIONS 

We have reviewed the transmission and localization 
wave properties of the complex disordered structures 
composed of (i) left-handed metamaterials, (ii) magneto- 
active optical materials, (iii) graphene superlattices, and 

(iv) nonlinear dielectric media. Interference origin of the 
wave localization, together with strong energy concen- 
tration, makes Anderson localization highly sensitive to 
weak modifications of the material properties. We have 
shown that the exotic properties of novel materials can 
drastically modify the main features of the wave localiza- 
tion. This brings about anomalous pronounced depen- 
dences of the wave transmittance and localization length 
on both wave and material parameters: frequency, angle 
of incidence, polarization, magnetization, nonlinearity, 
etc. As a result, remarkable phenomena appear, such as 
anti-(de-)localization, unidirectional transmission, slow- 
light propagation, and bistability. 

We have described a number of novel features ac- 
companying the wave localization in complex media, 
including: (i) dramatic suppression of localization in 
mixed stacks with left-handed metamaterials, (ii) Brew- 
ster, zero-£, and zero-/i delocalization, and (iii) anoma- 
lous transmission enhancement in periodic metamateri- 
als with only one disordered electromagnetic character- 
istics, (iv) nonreciprocal localization and unidirectional 
transmission through magneto-active disordered stacks, 

(v) angle-dependent transmission resonances in graphene 
superlattices, and (vi) bistability and nonreciprocity of 
transmission resonances in nonlinear disordered struc- 
tures. 

We believe that presented results significantly extend 
and enrich theory and potential application of the wave 
localization in complex disordered media. In particular, 
they provide a theoretical toolbox which can serve for de- 
sign of novel optical and electronic devices with unusual 
transport properties. 
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